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ABSTRACT 

We study the effects of turbulence on magnetic reconnection using three-dimensional direct numeri- 
cal simulations. This is the first attempt to test a model of fast magnetic reconnection in the presence 
of weak turbulence proposed by Lazarian & Vishniac (1999). This model predicts that weak turbu- 
lence, which is generically present in most of astrophysical systems, enhances the rate of reconnection 
by reducing the transverse scale for reconnection events and by allowing many independent flux recon- 
nection events to occur simultaneously. As a result the reconnection speed becomes independent of 
Ohmic resistivity and is determined by the magnetic field wandering induced by turbulence. We test 
the dependence of the reconnection speed on turbulent power, the energy injection scale and resistivity. 
We study the reconnection model with the open and experiment with the outflow boundary condi- 
tions and discuss the advantages and drawbacks of various setups. To test our results, we also perform 
simulations of turbulence with the same outflow boundaries but without a large scale field reversal, 
thus without large scale reconnection. To quantify the reconnection speed we use both an intuitive 
definition, i.e. the speed of the reconnected flux inflow, as well as a more sophisticated deflnition 
based on a formally derived analytical expression. Our results conflrm the predictions of the Lazarian 
& Vishniac model. In particular, we find that the reconnection speed is proportional to the square 
root of the injected power, as predicted by the model. The dependence on the injection scale for some 
of our models is a bit weaker than expected, i.e. compared to the predicted linear dependence 

on the injection scale, which may require some refinement of the model or may be due to the effects 
like finite size of the excitation region, which are not a part of the model. The reconnection speed was 
found to depend on the expected rate of magnetic field wandering and not on the magnitude of the 
guide field. In our models, we see no dependence on the guide field when its strength is comparable 
to the reconnected component. More importantly, while in the absence of turbulence we successfully 
reproduce the Sweet-Parker scaling of reconnection, in the presence of turbulence we do not observe 
any dependence on Ohmic resistivity, confirming that the reconnection of weakly stochastic field is 
fast. We also do not observe a dependence on anomalous resistivity, which suggests that the presence 
of anomalous effects, e.g. Hall MHD effects, may be irrelevant for astrophysical systems with weakly 
stochastic magnetic fields. 

Subject headings: galaxies: magnetic fields — physical processes: MHD — physical processes: turbu- 
lence — methods: numerical 



1. INTRODUCTION 

Magnetic fields play a key role in astrophysical pro- 
cesses such as star formation, the transport and accel- 
eration of cosm ic rays, a ccretion disks, solar phenom- 
ena, etc. (ICm tchcr 1991 lBeckll2002t ISchhckeiseill2004 
lElmegreen fc Scalo, .200 j ). Typically magnetic diffusion 
is very slow on astrophysical scales, so to a good ap- 
proximation we can treat magnetic fields as being purely 
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advected with the fiow, which is frequently referred to in 
the hterature as t he " frozen in" condition for the plasma 
(see lMoffatl[T978h . 

Do we expect the frozen in condition to be violated 
in typical astrophysical conditions? The answer to this 
question is a qualified yes. Indeed, most astrophysical 
flows are chaotic; when adjacent parcels of fluids do not 
move in the same direction, the frozen-in magnetic flelds 
become tangled. This can be easily visualized by view- 
ing magnetic flelds lines as threads moving with the fluid. 
When the distance between magnetic bundles of differ- 
ent direction becomes small, the finite resistivity of fluids 



2 



KOWAL ET AL. 



starts to be important. In a generic situation of 3D flows, 
the bundles of magnetic fields come into contact with 
their neighbors at an angle of the order unity. Over the 
small scales at which fluid resistivity is important, mag- 
netic field lines change their topology, or reconnect. How- 
ever, once magnetic field energies become large, bending 
magnetic fields on small scales requires energies much 
larger than the turbulent energies on that scale. The 
magnetic field lines become stiff, even in the presence of 
strong turbulence, and the scales associated with con- 
tact between regions with very different magnetic fields 
become large. 

What is the resulting reconnection speed? This 
is a vital question for many areas of astrophysics 
(see Biskamp '2000'; 'Priest & Forbedl2Q00|; iBhattachariei 
[2004: Zweib el k Yamada 2009 ). A first guess could be 
that magnetic reconnection is generically slow in astro- 
physical circumstances. There is a large disparity in the 
scales involved. The scale of the magnetic flux bundle 
is astronomically large. The microphysical scale over 
which the Ohmic dissipation is important is relatively 
small. In this case the crossing magnetic bundles will 
create ubiquitous unresolved intersections or "knots" in 
the fluid. Magnetic tension that arises from those inter- 
sections as the magnetic bundles press against each other 
should dramatically change the properties of magnetized 
fluids"'^. This would be devastating news for most nu- 
merical MHD simulations, as numerical diffusion in the 
simulations is high and magnetic bundles in simulations 
easily reconnect. 

To understand the difference between astrophysical re- 
connection and the one in numerical simulations, one 
should recall that the dimensionless combination that 
controls the reconnection rate is the Lundquist num- 
ber^, deflned as S' = LVa/ij, where L is the length of 
the reconnection layer (see Figure [U upper panel), Va is 
the Alfven velocity, and 77 is Ohmic diffusivity. Because 
of huge astrophysical sizes L involved, the astrophysi- 
cal Lundquist numbers are huge, e.g. for the ISM they 
are about 10^^, while present-day MHD simulations cor- 
respond to S" < 10*^. As the numerical efforts scale as 
L'^, where is the size of the box, it is not feasible at 
present and will not be feasible in the foreseeable future 
to have simulations with the sufficiently high Lundquist 
numbers. Incidentally, this also presents a problem for 
numerical simulations of magnetic reconnection unless 
one has theoretically-derived scaling relations to test. 
Even with the limited resolution, numerical simulations 
are a good tool to study scaling relations, the point that 
has been proved by successful numerical studies of MHD 
turbulence. 

Due to huge values of astrophysical Lundquist num- 
bers, any rate of reconnection that depends on S is ex- 
tremely slow^. Fast reconnection is reconnection that 
does not depend on resistivity. 

^ Don Cox (private communication) refers to this state of en- 
tangled magnetic fields as astrophysical Jello or astrophysical felt, 
to reflect the peculiar non-fluid properties of the hypothetical sub- 
stance. 

^ The magnetic Reynolds number, which is the ratio of the mag- 
netic field decay time to the eddy turnover time, is defined using 
the injection velocity vi as a characteristic speed instead of the 
Alfven speed Va, which is taken in the Lundquist number. 

^ The exception is the case when the dependence on S is loga- 
rithmic. 



What are the conditions that can make magnetic re- 
connection fast? Does it rely on special initial or bound- 
ary conditions for the flow, or does it require particular 
plasma effects? These are burning astrophysical ques- 
tions, which, for example, dcflne the extent that we 
can rely on numerical simulations of magnetized fluids 
as models of astrophysical phenomena. To understand 
processes of magnetic field generation associated with 
dynamos, the dynamics of the interstellar medium and 
accretion disks, or other related phenomena, it is im- 
portant to understand magnetic reconnection. In most 
cases, astrophysical reconnection is difficult to observe, 
with the notable ex ceptio n of solar flares (see ISturroca 
19661: iM asuda ct^aj] 1X9941 ) and gamma ray bursts (see 
Galama et al.lll998[ ). This sometimes creates an illusion 
that the importance of reconnection is limited to those 
selected phenomena. 

A famo us mode l of m agnet ic reconnectio n was sug- 
gested by iParked (|1957| ) and ISweetl (|1958l ). Unfortu- 
nately, this model, usually referred to as Sweet-Parker 
reconnection, provides very slow reconnection speeds. In 
this model the reconnection rate is inversely proportional 
to the square root of the Lundquist number. This very 
slow speed comes from a geometrical constraint. The 
current sheet dividing two magnetized regions must be 
very thin in order for Ohmic resistivity to be important, 
but the plasma trapped in the current sheet must follow 
the local magnetic field lines in order to escape. A narrow 
current sheet implies a highly restricted outflow. Natu- 
rally, if the Lundquist number is large, e.g. S = 10^^, 
the Sweet-Parker reconnection speed, Vsp ~ V^tS""^/^, 
is negligible. 

Fast reconnection has been investigated for many 
decades (seeEiskamp 2000l: iPriest fc Forbesl[2000l for re- 
views). In 1964 Petschek introduced the first fast mag- 
netic reconnection model (see .Petschek 19 64). He pro- 
posed that extended magnetic bundles come into contact 
over a tiny area determined by the Ohmic diffusivity. 
This configuration, called an X-point configuration, dif- 
fers dramatically from the expected generic configuration 
when magnetic bundles try to press their way through 
each other. Thus the first introduction of this model 
raised questions of dynamical self-consistency. An X- 
point configuration has to persist in the face of compres- 
sive bulk forces. However, numerical simulations have 
shown that an initial X-point configuration of magnetic 
field reconnection is unstable in the MHD limit for small 
values of the Ohmic diffusivity (|Biskamplll996f ) and the 
magnetic field will relax to a Sweet-Parker configuration. 
The physical explanation for this effect is simple. In the 
Petschek model shocks are required in order to maintain 
the geometry of the X-point. These shocks must persist 
and be supported by the flows driven by fast reconnec- 
tion. The simulations showed that the shocks fade away 
and the contact region spontaneously increases. 

More recently there has been progress in understand- 
ing the role of coUisionless modes in stabilizing X-point 
reconnection. In particular, simulations usin g Hall MHD 
have produ ced stable X-point reconnection (jShav et al.l 
fT99i. I200I . This was the flrst numerical demonstra- 
tion of fast reconnection, which raised the hope in astro- 
physical community that the problem of astrophysical 
reconnection can be solved. There is continuing discus- 
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sion of whether this kind of fast reconnection persists at 
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12004') yielding reconnection rates that are not fast, but 
depend on resistivity. More importantly, one may also 
wonder to what extent Petscheck-type collisionless recon- 
nection solves the problem of describing magnetic field 
dynamics in astrophysical settings. Apart from the issue 
of how natural it is to produce conditions that lead to 
Petscheck reconnection, one should note that the require- 
ment that the upper limit on the collision rate required 
by this proc ess is ex t remel y restrictive. For instance, 
estimates in lYamadal (|2007t ) show that the scale of the 
reconnection current sheet should not exceed approxi- 
mately 40 electron mean free paths. This is not satisfied 
in many astrop hysical environ ments including the inter- 
stellar medium (|Yamadall2007D . 

What happens in a collisional medium? Is reconnec- 
tion in the interstellar medium slow? The latter, as men- 
tioned above, would have catastrophic implications for 
the entire current crop of MHD simulations of interstel- 
lar processes. 

Nearly simultaneously with the discovery of collision- 
less X-point reconnection, a model of magnetic recon- 
nection in the presen ce of a weakly stochastic m agnetic 
field was proposed bv lLazarian fc Vishniad (|1999l hence- 
forth LV99). They claimed that the laminar magnetic 
fields considered in both Petscheck and Sweet-Parker re- 
connection models are exceptional in astrophysics. If we 
consider the interstellar medium as an example, a so- 
called " Big Power Law in the Sky" indicates the presence 
of turbulence o n scales from tens of pa rsecs to thousands 
of kilometers (jArmstrong et al.lll995[ ). Thus, the fact 
that reconnection may be slow in a collisional fluid with 
laminar field lines does not necessarily imply that recon- 
nection is generically slow. The model in LV99 predicted 
extended current sheets and wide outflows limited by the 
diffusive spread of magnetic field lines. In the limit of 
very weak turbulence or very large resistivities the model 
converges to the Sweet-Parker model. 

For many years the quantitative predictions of the 
stochastic reconnection model in LV99 have not been 
tested explicitly. There is some implicit evidence in the 
favor of this model, e.g. observations of the thick re- 
connection cuTrent_outflOT;_j^£i^ in the So- 
lar flares ([Ciaravella fc Ravmondl [2008^ . Such outflow 
regions are incompatible with Petscheck reconnection. 
This is suggestive, but falls far short of constituting a 
definitive proof of the LV99 model. 

It is, in general, extremely difficult to us e analytic tools 
to tes t turbulent processes. However, lEvink fc Aluid 
(|2006l ) have studied the necessary conditions for recon- 
nection in a turbulent conducting medium in the limit of 
infinitesimal resistivity. They list three separate neces- 
sary conditions, and prove that the flux conservation can 
be violated at an instant of time for an arbitrarily small 
scale, if at least one of the conditions is satisfied. The one 
that concerns us here is that there must be overlapping 
vortex and current sheets in the fluid. This is trivially 
satisfied if turbulent eddies on all scales give rise to cur- 
rent sheets with inflow/outflow properties of the kind we 
study here. While this shows that the LV99 model is 
consistent with the physical requirements for fast recon- 



nection, it does not prove that the model is viable. 

Numerical simulations have provided the main test- 
ing ground for collisionless X-point reconnection. These 
simulations are typically performed in 2D with peri- 
odic boundary conditions and cannot be continued for 
longer than one dynamical time. This approach is in- 
appropriate for stochastic reconnection where the goal is 
to obtain a stationary reconnection rate averaged over 
many dynamica l times . For this reason we have followed 
iDaughton eFall (|2006[ ) in using open outflow boundary 
conditions. Moreover, the LV99 model of reconnection 
is intrinsically three dimensional. A simplification, how- 
ever, arises from the fact that LV99 predicts that fast re- 
connection takes place in the MHD approximation, mak- 
ing plasma effects irrelevant. 

In this paper we describe the choice of boundary condi- 
tions as well as an appropriate measure of reconnection. 
Indeed, while the rate of reconnection can be trivially 
measured in the case of laminar magnetic fields, addi- 
tional care is required in the estimation of the reconnec- 
tion rate in the presence of a stochastic magnetic field. 

In fj2] we review the LV99 model of reconnection and 
its theoretical predictions. In Sj3]we describe in detail the 
numerical model we have built and studied in this paper. 
In 21 '^f' present methods used throughout the paper to 
measure the reconnection rate and we introduce a new, 
more general, method of estimating the reconnection rate 
which includes all processes contributing to the change of 
magnetic flux. In fj5]we present an extensive description 
of results obtained from studying our numeric model, 
which we discuss later in fJS) In fJ7]we set forth our main 
conclusions. 

2. LAZARIAN-VISHNIAC (1999) MODEL 

We begin by reviewing the b asic features of the 
stocha stic reconnection model in iLazarian fc Vishniad 
(|1999| ). This model might be seen as a generalization of 
the Sweet-Parker model (see Fig.[T|) for the case of turbu- 
lent fields. Like the Sweet-Parker model it deals with an 
extremely generic configuration, which should arise nat- 
urally as magnetic flux tubes try to make their way one 
through another. This avoids the problems related to the 
preservation of outflow conditions which plague attempts 
to explain magnetic reconnection via Petscheck-type so- 
lutions. For example, if the outflow of reconnected flux 
and entrained matter is temporarily impeded through a 
fluctuation in the outflow, reconnection will slow down, 
but will not permanently change its nature. 

The essential difference between the Sweet-Parker 
model and the model in LV99 is that while in the Sweet- 
Parker model the outflow is limited by microphysical 
Ohmic diffusivity, in the LV99 model the large-scale mag- 
netic fleld wandering determines the thickness of outflow. 
For extremely weak turbulence, when the range of mag- 
netic fleld wandering becomes smaller than the width of 
the Sweet-Parker layer LS~^^^, the reconnection rate re- 
duces to the Sweet-Parker rate^. 

* For the sake of simplicity, we do not raise here issues related to 
the stability of the Sweet-Parker current sheet. Estimates in LV99 
suggest that laminar current sheets subject to tearing instabilities 
have reconnection rates which are just a bit faster than in the 
original Sweet-Parker model. These instabilities may, however, be 
important for increasing the 3D stochasticity of magnetic field lines 
and thus initiating fast reconnection. In this paper we do not 
address this effect either. 
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Fig. 1. — Upper plot: Sweet-Parker model of reconnection. The 
outflow is limited to a thin width 5, which is determined by Ohmie 
diffusivity. The other scale is an astrophysical scale L ^ 5. Mag- 
netic field lines are assumed to be laminar. Middle plot: Turbulent 
reconnection model that accounts for the stochasticity of magnetic 
field lines. The stochasticity introduced by turbulence is weak and 
the direction of the mean field is clearly defined. The outfiow is 
limited by the diffusion of magnetic field lines, which depends on 
macroscopic field line wandering rather than on microscales deter- 
mined by resistivity. Low plot: An individual small scale reconnec- 
tion region. The reconnection over small patches of magnetic field 
determines the local reconnection rate. The global reconnection 
rate is substantially larger as many independent patches reconnect 
simultaneously. Conservatively, the LV99 model assumes that the 
small scale events ha ppen at a slow Sweet-Parker rate. Following 
ILazarian eFaTI I I2004I) . 

LV99 consider the case of a large scale, well-ordered 
magnetic field, of the kind that is normally used as a 
starting point for discussions of reconnection. In the 
presence of turbulence we anticipate that the field will 
have some small scale 'wandering'. On any given scale the 
typical angle by which field lines differ from their neigh- 
bors is (j> <^ I, and this angle persists for a distance along 
the field lines Ay with a correlation distance Aj_ across 
field lines (see Fig. [T]) . 

LV99 suggested that the presence of a random mag- 
netic field component substantially enhances the recon- 
nection rate, leading to fast reconnection. There are two 
phenomena mainly responsible for this: 

• only a small fraction of any magnetic field line is 
subject to direct Ohmic annihilation. The fraction 
of magnetic energy that goes directly into heating 
the fluid drops down to zero as the fluid resistivity 
vanishes, and 

• the presence of turbulence enables many magnetic 
field lines to enter the reconnection zone simulta- 
neously. 

In order to do quantitative estimates one has to 
adopt a model of MHD turbulence. Attemp t to con- 
struct such a model ca n be tr aced back to ilroshnikovl 
(1961 and iKraichnanI (fT965h papers as well as to 



later fundamental works (e.g. [Montgomery &: Turneil 
fl98lt IShebalin et al.l[l98l IHigdod HQS'T see also book 
by Biskamp 2003). A model for realistic compress- 
ible fluids may be constr ucted on the basis of the 
(jGoldreich fc Sridhad [l995l . henceforth GS95) model of 
incompressible turbulence^. Reasonable alternatives are 
expected to have only marginal impact on the nature 
of field line stochasticity (see als olCho Sz V ishniac"2000j: 
Maron fc G oldrcich 2001; Lith wick fc Go ldrcich 20011; 
Cho et al.l 120021 : iCho fc Laza rian 2001 l2003ll . provided 



that the GS95 scalings are generalized to include the 
case of weak MHD turbulence. The equations in LV99 
reflect the fact that when the turbulence is injected at 
subAlfvenic velocities, the cascading is initially weak and 
the perturbations can be reasonably well represented by 
a collection of Alfven waves undergoing occasional inter- 
actions. However, the strength of the non-linear interac- 
tions increases as the perpendicular scale decreases and 
at the scale l± = 1{vi/VaY the cascade becomes strong, 
i.e. with the GS95 critical balance between motions at 
parallel and perpendicular scales l±/vj_ ^ 1/Va satisfied. 
The corresponding velocity at l± is v± ^ V^a(i'//Va)^ 
where vi is the velocity at which strong turbulence is 
being injected. For vi <^ Va this scale may be small 
and the field wandering induced by the turbulence is re- 
duced^. The characteristics of the strong turbulent cas- 
cade dominate transport phenomena and predictions of 
reconnection speeds. In contrast the weak turbulent cas- 
cade has very little effect on field line stochasticity, since 
in this regime motions are largely periodic. 

The modification of the global constraint induced by 
mass conservation in the presence of a stochastic mag- 
netic field component is self-evident. Instead of being 
squeezed from a layer whose width is determined by 
Ohmic diffusion, the plasma diffuses through a much 
broader layer, S ~ (y^)^^^ (see Fig. [1]), determined by 
the diffusion of magnetic field lines. This suggests an up- 
per hmit on the reconnection speed of ~ Va 
This will be the actual speed of reconnection if the 
progress of reconnection in the current sheet does not 
impose a smaller limit. The value of (y^)^^^ can be deter- 
mined once a particular model of turbulence is adopted, 
but it is obvious from the very beginning that this value 
is determined by field wandering rather than Ohmic dif- 
fusion. 

Following ILazarian et al.l (|2004 hereafter LVC04) we 
can generalize this upper limit to include the whole range 
of scales between the large scale eddy size and the current 
sheet thickness. Consider two points initially separated 

^ The exact scalings of the MHD turbulence are still subject 
to debate and different corrections to the original GS95 scalings 
have been proposed to account for dynamical alignment, polariza- 
tion intermittenc y and non-l ocality of the cascade (see ]3oldvre3 
2005, 2006; Bercs nvak fc Laza rian 2006, 2008b; Gogoberidze 2007:1). 
However, calculations in LV99 show that reconnection rates are 
only weakly dependent on the exact model of turbulence. Subtle 
refinements of the scalings that are currently debated can only be 
of marginal significance for this paper. For the sake of simplic- 
ity, we adopt in what follows the original GS95 scaling. We also do 
not cons ider imbalanccd MHD turbulence (sec Lithwick et al. 200^; 
IChandran.i2008. : .Beresnvak fc Lazariann2008a . .200 9). although for 
the modest degrees of imbalance between the turbulent energy 
fluxes moving in the opposite directions, we do not expect sub- 
stantial changes in our results 

^ These comments ignore the compressible modes, but their role 
in field wandering is marginal anyhow (LV99). 
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by the thickness of the current sheet. Translating both 
points along the field lines we can define their rms separa- 
tion as a function of distance along the field lines, (5(A||). 
This gives an upper limit on the local reconnection speed 
which is just ^ {5/\\\)Va- However, the large scale cur- 
rent sheet contains many such reconnection regions, each 
involving the reconnection of independent field lines. In 
order to estimate the global reconnection rate we need 
to multiply each local reconnection rate with the number 
of simultaneous reconnection events happening on that 
scale. Taking L/\\\ as our best estimate we get a global 
limit of 

-5(A||) L 



Vrpr < niin 



All Ai 



■Va 



(1) 



which should be evaluated for all A|| between the length 
of the current sheet and the length of an individual piece 
of the current sheet such that translation along that 
piece will increase the rms separation by the thickness 
of the current sheet. As long as we are in the regime 
of strong turbulence this expression has a rather sim- 
ple form. Each eddy scatters field lines by roughly its 
own width, so an eddy of parallel length A|| will have a 
corresponding 6 ^ A_l(A||). With this in mind Eq. ((1]) 
becomes 



Vrpr rnin 



-^Va 



(2) 



At its minimum Eq. ^ for GS95 model provides Vrec < 
Va, which is a natural limit of the reconnection speed 
and not a constraint. In the Goldreich-Sridhar model 

2/3 

of turbulence A|| oc A_l , so the minimum value of this 
expression corresponds to the largest scales. (Once A|| 
exceeds the scale of the largest eddies we change models 
to a random walk with a fixed step size.) 

This should make it clear why we are insensitive to 
the exact model of strong turbulence used in our calcu- 
lation. In order to get a different answer we would need 

1/2 

A|| cx A_L or less, which is not supported by any of the 
available numerical evidence. On the other hand, how 
one models the weak turbulence does have an impact on 
the reconnection rate, since it will define the length of 
the largest scale eddies. Also, the presence of strong vis- 
cous damping can complicate our argument, by creating 
a laminar field line regime on very small scales and reduc- 
ing S far below its value in the strong turbulent regime. 
This point is discussed at length in LVC04 where it was 
argued that in some of the colder and denser phases of 
the interstellar medium we should expect a reduction in 
the reconnection speed by an order of magnitude or more. 
The numerical work presented here includes only a min- 
imal amount of viscosity, so testing these arguments is 
beyond the scope of this paper. 

The basic hypothesis of LV99 is that the the speed of 
magnetic reconnection in the presence of 3D stochastic- 
ity of magnetic field lines is given by the most stringent 
constraint imposed by the bottleneck condition expressed 
in Eq. (P) With the GS95 model of turbulence LV99 

''' In LV99 other processes that can impede reconnection were 
found to be less restrictive. For instance, the tangle of reconnec- 
tion field lines crossing the current sheet will need to reconnect 
repeatedly before individual flux elements can leave the current 
sheet behind. The rate at which this occurs can be estimated by 



obtained: 



Va min 



1/2 



1/2- 



where / and vi are the energy injection scale and tur- 
bulent velocity at this scale respectively. Note that the 
combination Va{vi/Va)'^ is the velocity Vtrans, i-e. the 
velocity at which the cascade transfers to the strong 
regime. This refiects the importance of strong turbulence 
for magnetic field wandering and magnetic reconnection. 
The term (L//)^/^ reflects the random walk which takes 
place for L < I Evidently, the goal of our paper is to test 
formula ([3]). 

There are several important characteristic features of 
Eq. ([3]). First, and most important, there is no depen- 
dence on resistivity. Second, in general we expect re- 
connection to be fast since in most cases the parameter 
ratios that enter the expression, i.e. the length of the 
reconnection layer L divided by the injection scale /, and 
the injection velocity vi divided by the Alfven velocity 
Va are of order unity. Finally, we note that in partic- 
ular situations when turbulence is extremely weak the 
reconnection speed may be small. 

Given the limited dynamical range of the simulations, 
we are forced to inject turbulent energy on scales less 
than L. Also, it is easier to control not vi, but the energy 
injection power P. The power in the turbulent cascade 
is P ^ Vf^i^^ij{VA / 1) or vf/{lVA), which reflects the fact 
that the turbulence is weak at the injection scale I (see 
LV99). Using this in Eq. §^ we get 



Vr, 



(4) 



which is the prediction we will test here. In what follows 
we refer to the injection power and scale using Pj„j and 
Unj, respectively. 

3. NUMERICAL SETUP 

3.1. Governing Equations 

We use a higher-order shock-capturing Godunov- 
type scheme based on the essentially non oscilla- 
tory (ENO) spacial reco nstruction and Runge-Kutta 
(RK) time integration (see iLondrillo fc Del Zaniiall200Cll : 
iDel Zanna et a]|f2003L e.g.) to solve isothermal non-ideal 
MHD equations, 

|^ + V.(pi;)=.0, (5) 



dpv 



V • 



pvv - 



El 



/ BB 

4tt 



= /■ (6) 



f + E^O. (7) 

where p and v are plasma density and velocity, respec- 
tively, A is the vector potential, E = —vx B + rjj is the 
electric field, B=VxAis the magnetic field, j = V x B 
is the current density, a is the isothermal speed of sound, 
rj is the resistivity coefficient, and / represe nts the forcing 
term. We used Harten-Lax-van Leer (HLL. lHarten et aLl 

assuming that it constitutes the real bottleneck in reconnection 
events, and then analyzing each flux element reconnection as part 
of a self-similar system of such events. This turns out not to impede 
reconnection. 
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XZ plo 
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Fig. 2. — Schematic 3D visualization of the reconnection problem 
setup. Upper and lower parts of the domain contain incoming 
oppositely directed magnetic field lines (black and grey solid lines, 
respectively). The inflow bends the lines toward the center where 
they enter the diffusion region. In the diffusion region the lines 
reconnect and their product (dashed lines) is ejected along the X 
direction. 

119831 ) and HLLD Riemann solvers (jMignond [20071 ) for 
solving the isothermal MHD equations. These schemes 
have widely different dissipation properties. The HLL 
solver was initially developed for solving Euler equations 
and later adopted for solution of general set of time de- 
pendent differential equations by averaging the Riemann 
fan over a region bound by the minimum and maximum 
local characteristic speeds. In this way one intermedi- 
ate state was constructed. The HLLD Riemann solver 
takes into account the discontinuities resulting from the 
presence of magnetic field as well, separating this inter- 
mediate state into multiple intermediate states, resolving 
e.g. Alfven waves with much less dissipation. In this pa- 
per we are considering the quasi-incompressible regime, 
where most of energy is transported by Alfven waves, 
therefore the application of the HLLD solver seems to be 
a better choice. We incorporated the field interpolated 
constrain ed transport (CT) s cheme based on a staggered 
mesh (see iLondrillo fc Del Z anna 2000) into the integra- 
tion of the induction equation (Eq. [7]) to maintain the 
V • B = constraint numerically. 

Some selected simulations that we perform include 
anomalous resistivity modeled as 



ri = ilu+ Va 



-I H 



Jcrit 



Jcrit 



(8) 



where r]u and rja describe uniform and anomalous resis- 
tivity coefficients, respectively, jcrit is the critical level 
of the absolute value of current density above which the 
anomalous effects start to work, and H isa. step function. 
However, for most of our simulations rja = 0. 

3.2. Model Description and Initial Conditions 

Figures [2] and [3] show a 3D visualization and 2D pro- 
jections of the reconnection problem setup. The domain 
contains two regions of oppositely directed magnetic field 
lines (see Fig. [2] and the left panel of Fig. [3]). The incom- 
ing lines (solid lines in Fig. [2|) are bent by the inflow Vi„ 
and enter the diffusion region. The diffusion region is 
characterized by the longitudinal scale A and its thick- 
ness 6 (see the left panel of Fig. [3|) . The diffusion region 
extends along the full Z direction of the system. The in- 
coming magnetic lines are not perfectly antiparallel. The 



projection of the magnetic topology on the XZ plane 
shows that the lines in the upper region (solid lines in 
the right panel of Fig. [3]) and in the lower region (dashed 
lines) create an angle a determined by the strength of 
the shared component Bqz . Once the incoming magnetic 
lines enter the diffusion region, they are reconnected and 
the product of this process is ejected along X direction 
with a speed Vout (see Fig.[2]and the left panel of Fig. [3]). 

Our initial magnetic field is a Harris current sheet of 
the form Bx{x, y, z) = Bq^ tanh(?//0) initialized using the 
magnetic vector potential Az{x,y,z) = In | cosh(?//6')|. 
In addition, we use a uniform guide field Bz{x,y, z) = 
Bqz = const. The initial setup is completed by setting 
the density profile from the condition of the uniform to- 
tal (thermal plus magnetic) pressure prit = 0, x, y, z) = 
const and setting the initial velocity to zero everywhere. 

In order to initiate magnetic reconnection we 
add a small initial perturbation of vector potential 
5Az{x,yTz) ~ (5i?0K cos(27ra;) exp[— (y/d)^] to the initial 
configuration of Az{t = 0,x,y, z). Parameters (5i?ox and 
d describe the strength of the initial perturbation and 
thickness of the perturbed region, respectively. 

We use dimensionless equations, so that the strength 
of the magnetic field is expressed in terms of the Alfven 
velocity (defined by the antiparallel component of mag- 
netic field) and the unperturbed density = 1. All 
other velocities are expressed as fractions of the fiducial 
Alfven speed. The length of the box in X direction de- 
fines the unit of distance and time is measured in units 
of Lx/Va- Initially, we set the strength of antiparallel 
magnetic field component Bq^ to 1.0 and we vary the 
guide field Bqz between 0.0 and 1.0, which corresponds 
to the range of angle a G (0°, 90°). For a particular run 
we use Bqz = 4.0, which corresponds to a ~ 152°. The 
speed of sound is set to 4.0 to suppress compressibility 
in the system. In order to study the importance of resis- 
tivity in the reconnection process we vary the resistivity 
coefficient rju between values 5 • 10~^ and 5 • 10~^ which 
are expressed in the dimensionless units. In the models 
where we include anomalous effects, we vary the anoma- 
lous resistivity coefficient r]a between 0.0 and 2 • 10^''. 
The parameters describing the initial perturbation are 
set to SBqx = 0.05 and d = 0.1. 

3.3. Boundary Conditions 

Our numerical model of LV99 reconnection evolves in 
a box with dimensions = Lz = 1 and Ly = 2 with 
resolution 256x512x256. It is extended in the Y direction 
in order to move the inflow boundaries far from the injec- 
tion region. This minimizes the influence of the injected 
turbulence on the inflow. 

We use three different types of boundary conditions 
in our models, each of them set along different direc- 
tion: outflow boundary conditions along X direction, in- 
flow boundary conditions along Y direction and periodic 
boundary conditions along Z direction. 

We tested different types of open boundary conditions, 
including perfectly permeable wave-type boundary con- 
ditions derived from the linearized M HD equations sim - 
ilar to the method of characteristics (H edstromI Il979l ). 
which seem to be the right choice. However, they are 
expensive since they require to track the linear MHD 
waves crossing the boundaries and they have a relatively 
small advantage over the simple "zero-gradient" bound- 
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open or periodic 




open or periodic 



Fig. 3. — 3D magnetic field configuration in the studied problem projected on the XY (left) and XZ (right) planes. Left: XY projection 
of the magnetic field lines. The grey area describes the diffusion region where the incoming field lines reconnect. The longitudinal and 
transverse scales of the diffusion region are described by the parameters A and S, respectively. represents the direction of the guide field 
with respect to the plane of projection. We use inflow and outflow boundary conditions at X and Y directions, respectively. Right: XZ 
projection of the magnetic field lines as seen from the top. Solid and dashed lines show the incoming field lines from the upper and lower 
parts of the domain, respectively. We see that the oppositely directed field lines are not perfectly antiparallel but they create an angle a 
determined by the strength of the shared component B^. The Z boundary conditions are open or periodic, depending on the model. 



ary conditions in the application to our models, which we 
found the simplest and most robust for our purposes. In 
this type of open boundary conditions we set the normal 
derivatives of the fluid variables (density and momen- 
tum) to zero. This guarantees that all waves generated 
in the system are free to leave the box without significant 
reflections. In the presence of turbulence, however, this 
type of boundary does not enforce a constant mean value 
of density at the inflow boundary, which can result in a 
small loss of the total mass in the system. We will show 
later that the total mass loss is very small and does not 
influence our results. 

As for the vector potential, we have tested 
Daughton's approach for open boundary conditions (see 
iDaughton et al.l l2006f) in which the authors solved the 
properly set wave equations for each component. How- 
ever, this approach was limited to the 2.5 dimensional 
case and its extension to the full 3D case requires a sig- 
nificant increase of complexity and, in our opinion, re- 
sults in an unacceptably ex pensive method. Moreover, 
since IDaughton et al.l (|2006f ) studied the problem of re- 
connection using a kinetic approach, they used the speed 
of light to evolve traveling electromagnetic waves. In 
the MHD approximation, this would imply a very large 
reduction in the time step, since the stability of the 
scheme de pends on the Courant-Friedrichs-Lewy (CFL, 
ICouranFT t al. 1928) condition, which in turn depends on 
the maximum speed in the system. This issue excluded 
this kind of approach from our consideration. 

In the end we decided to use the simplest approach. 
We evolve the induction equation for the vector poten- 
tial A based on a staggered mesh, which means that the 
vector potential components are located in the edge cen- 
ters, i.e. A,4i,j + 1/2, k+ 1/2), Ay{i + l/2,j,k + l/2), 
and Az{i + 1/2, j + 1/2, k), where {i,j,k) describes the 
ce ll center position. Thus , using the method described 
in iLondrillo &: Del Zannal ()2000D we calculate magnetic 
field components located at the interface centers, [i + 
1/2, J, k), By(i,j + l/2,k), and B,(i, j, /c + 1/2). These 
locations of collocation points for the vector potential 
completed by a CT integration of the induction equa- 
tion guarantees the divergence- free evolution of magnetic 
field. 



In the treatment of vector potential A at the bound- 
ary we set its components transverse to the considered 
boundary using the first order extrapolation, while the 
normal derivative of the normal component is set to 
zero. In this way the normal derivatives of transverse 
components of magnetic field are zero, while the normal 
component of magnetic field is calculated from the zero- 
divergence condition V ■ B — Q. 

This approach prevents from the generation of diver- 
gence of B at the boundaries, however, it creates a small 
jump of the fluxes in the momentum equation across 
the boundary resulting from the presence of non-zero 
terms {—Bx,By,Bz)dxBx at the X outflow boundary 
and {Bx, —By, Bz) dyBy at the Y inflow boundary. In 
order to justify the importance of these terms we have 
estimated directly from our models the velocity incre- 
ment they produce at each time step. In models with the 
strongest turbulence the magnitudes of these terms were 
of order of 10"^ and 10"^ of Va at the X and Y bound- 
aries, respectively, signifying their negligible importance 
in the presence of strong inflow and outflow which are of 
order of Alfven speed in our models. 

There is another constrain on the boundary conditions 
applying to the induction equations resulting from its 
non-ideality, i.e. the presence of resistive term and the 
treatment of current density j at the boundaries. When 
we set normal derivatives of transverse components of 
magnetic field to zero, we violate the continuity of cur- 
rent density. In the ideal case, ry = 0, the resistivity term 
is unimportant even if the current density is not contin- 
uous. However, once we introduce a non-zero resistivity 
coefficient, the resistive term rjjm the induction equation 
becomes discontinuous at the boundary and we introduce 
a jump of electromotive force £ = v x B — rjj across the 
boundary. For example, considering the Z component of 
electromotive force = v^By — VyB^ — rj {dxBy — dyB^) 
we obtain a jump of across the X boundary equal to 
^x£z = {£1 — £z) = V^xBy, where indices d (domain) 
and b (boundary) stand for the limit values from the 
left and right sides of the considered boundary, respec- 
tively. This results in a constant generation of By at the 
boundary if the second derivative d^By is not zero at 
that location. One way to remove this undesired effect is 
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by setting the second order normal derivative of the field 
components to zero. However, there is no straightforward 
way to do this since it introduces another constraint on 
the boundary values and requires complex and expen- 
sive methods. We use another way of diminishing the 
importance of resistive terms and current density jump 
at the boundary. Instead of modifying the boundary con- 
ditions, we introduce a zone of decaying resistivity close 
to the boundary. This means that the resistivity is con- 
stant in most of the computational domain. Only in a 
thin zone near the boundary, the value of resistivity rju 
decays down to a very small value, which is estimated to 
be around the numerical resistivity ?7„ of our code. In 
our models we adopt the value of r]n ~ 3 ■ 10^"*. Nat- 
urally, by introducing such resistivity we do not change 
the reconnection speeds. There are two advantages to 
this approach. First, the outflow in our models is not 
constrained by the jump of resistive terms growing with 
the increasing value of resistivity and is determined 
solely by the evolution in the central part of the domain, 
where the resistivity is constant and uniform. Second, we 
control the importance of numerical effects by introduc- 
ing a small, but larger than numerical resistivity value 
of ry„ at the boundary. The validation of this method is 
presented in ^ 15.2.41 

The above boundary conditions represent open bound- 
aries, which adjust during the evolution of the system. 
This means that we do not set fixed values of the fluid 
and magnetic variables and do not drive the flow at the 
boundaries in order to achieve a stationary reconnection 
process in the system. 

3.4. Model of Turbulence 

In our m odels we drive t urbulence using a method de- 
scribed bv lAlveliui (|1999( ). The forcing is implemented 
in spectral space where it is concentrated around a wave 
vector kinj corresponding to the injection scale linj. We 
perturb a number Nf of discrete Fourier components 
of velocity in a shell extending from kinj — ^kinj to 
kinj + Afcinj with a Gaussian profile of the half width 
kc and the peak amplitude Vf at the injection scale (see 
Tab. [1] for exact values in all models). In all models 
kc = 0.4. Since we can control the scale of injection, the 
power input is introduced into the flow at an arbitrary 
scale. The amplitude of driving is solely determined by 
its power Pi„j , the number of driven Fourier components 
and the time step of driving Atf, which in2 all our mod- 
els is equal to 10~^. The parameters describing our forc- 
ing do not change during the evolution of the system. 
Because we perturb discrete number Nf of Fourier com- 
ponents which depends on the injection scale kinj and 
the thickness of perturbed shell Akinj , the amplitude v f 
varies with the injection scale in order to keep the same 
power input Pinj- In our models w/ is always a small 
fraction of Alfven speed (see Tab.[T]). 

The randomness in time makes the force neutral in 
the sense that it does not directly correlate with any of 
the time scales of the turbulent flow, and it also deter- 
mines the power input solely by the force-force correla- 
tion. This means that it is possible to generate different 
desirable turbulence states, such as axisymmetric turbu- 
lence, where the degree of anisotropy of the forcing can 
be chosen a priori through the forcing parameters. In the 
models presented in this paper we use isotropic forcing 



only 

In particular, the total amount of power input from 
the forcing can be set to balance a desired dissipation at 
a statistically stationary state. In order to contribute to 
the input power in the discrete equations from the force- 
force correlation only, the force is determined so that the 
velocity-force correlation vanishes for each Fourier mode. 
The procedure of reduc i ng th e velocity-force correlation 
is described in lAlveliud (| 19990 . 

We drive turbulence in a subvolume of the domain. 
The size of the subvolume is determined by two scales, 
the radius on the XZ plane around the center of the 
domain and the height hd describing the thickness of the 
driving region from the midplane. In this way we avoid 
driving turbulence at the boundary and reduce the influ- 
ence of driving on the inflow or outflow. All models are 
evolved without turbulence for several dynamical times 
in order to allow the system to achieve laminar stationary 
reconnection. Then, at a given time th we start driving 
turbulence, increasing its amplitude to the desired level, 
until time te- In this way we let the system to adjust to 
a new state. From time ie the turbulence is driven with 
the full power Pinj ■ 

On the right hand side of Eq. ^ , the forcing is repre- 
sented by a function f— pa, where p is local density and 
a is random acceleration calculated using the method 
described above. 

We do not set the viscosity coefficient explicitly in our 
models. The scale at which the dissipation starts to be 
important is defined by the numerical diffusivity of the 
scheme. The ENO -type schemes are considered relatively 
low diffusion (see iLiu fc OsheTI Il998t iLevv et all Il999i . 
e.g.). The numerical diffusion depends not only on the 
adopted numerical scheme but also on the "smoothness" 
of the solution, so it changes locally in the system. In 
addition, it is also a time-varying quantity. All these 
problems make its estimation difficult and incompara- 
ble between different applications. However, the dissi- 
pation scales can be estimated approximately from the 
velo city spectra. Supported by our studies of tur bulence 
(see iKowal et all [2007; Kowal & Lazariaii I2007D where 
we used similar code, we estimated the dissipation scale 
k^ « 30 for the resolution of models presented here. 

3.5. Table of Simulated Models 

In Table [1] we list parameters of all the models pre- 
sented in this paper. We divided them into several 
groups. In each group we calculated models in order to 
study the dependence of the reconnection rate on a char- 
acteristic parameter of turbulence or resistivity. We have 
studied the dependence of reconnection on the power 
of turbulence (models "PD"), injection scale (models 
"SD"), uniform resistivity (models "RD"), anomalous re- 
sistivity (models "AD"), and dependence on the guide 
field (models "BD"). 

Among all parameters of the model we list those which 
vary, i.e. the strength of guide field i?oz, the uniform 
and anomalous resistivities, rju and rja, respectively, the 
power of turbulence Pinj and its injection scale kinj with 
the half-thickness of the injection shell Akinj, the num- 
ber of perturbed Fourier components of velocity Nf and 
the amplitude of perturbation Vf at the injection scale. 
In addition, we include some of the most important 
parameters obtained from the simulations, such as the 
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TABLE 1 
List of models. 



Name 




Vu [10-3] 


Va [10-3] 


Pinj 






Nf 


Vf 


vi (t=12) 


Va 


/3 = p/pmag 


PD 


0.1 


1.0 


0.0 


0.1 


8 


0.5 


96 


0.0015 


0.031 


1.05 


31.7 




0.1 


1.0 


0.0 


0.2 


8 


0.5 


96 


0.0022 


0.041 


1.05 


31.7 




0.1 


1.0 


0.0 


0.5 


8 


0.5 


96 


0.0035 


0.051 


1.05 


31.7 




0.1 


1.0 


0.0 


1.0 


8 


0.5 


96 


0.0049 


0.065 


1.05 


31.7 




0.1 


1.0 


0.0 


2.0 


8 


0.5 


96 


0.0069 


0.084 


1.05 


31.7 




1.0 


1.0 


0.0 


0.1 


8 


0.5 


96 


0.0015 


0.042 


1.41 


16.0 




1.0 


1.0 


0.0 


0.2 


8 


0.5 


96 


0.0022 


0.043 


1.41 


16.0 




1.0 


1.0 


0.0 


0.5 


8 


0.5 


96 


0.0035 


0.056 


1.41 


16.0 




1.0 


1.0 


0.0 


1.0 


8 


0.5 


96 


0.0049 


0.071 


1.41 


16.0 




1.0 


1.0 


0.0 


2.0 


8 


0.5 


96 


0.0069 


0.083 


1.41 


16.0 


SD 


0.1 


1.0 


0.0 


1.0 


5 


0.5 


24 


0.0028 


0.092 


1.05 


31.7 




0.1 


1.0 


0.0 


1.0 


8 


0.5 


96 


0.0049 


0.063 


1.05 


31.7 




0.1 


1.0 


0.0 


1.0 


12 


1.0 


48 


0.0018 


0.059 


1.05 


31.7 




0.1 


1.0 


0.0 


1.0 


16 


1.0 


96 


0.0049 


0.035 


1.05 


31.7 




0.1 


1.0 


0.0 


1.0 


25 


1.0 


265 


0.0014 


0.030 


1.05 


31.7 




1.0 


1.0 


0.0 


1.0 


5 


0.5 


24 


0.0028 


0.092 


1.41 


16.0 




1.0 


1.0 


0.0 


1.0 


8 


0.5 


96 


0.0049 


0.068 


1.41 


16.0 




1.0 


1.0 


0.0 


1.0 


12 


1.0 


48 


0.0018 


0.059 


1.41 


16.0 




1.0 


1.0 


0.0 


1.0 


16 


1.0 


96 


0.0049 


0.034 


1.41 


16.0 




1.0 


1.0 


0.0 


1.0 


25 


1.0 


265 


0.0014 


0.027 


1.41 


16.0 


RD 


0.1 


0.5 


0.0 


1.0 


8 


0.5 


96 


0.0049 


0.066 


1.05 


31.7 




0.1 


1.0 


0.0 


1.0 


8 


0.5 


96 


0.0049 


0.065 


1.05 


31.7 




0.1 


2.0 


0.0 


1.0 


8 


0.5 


96 


0.0049 


0.067 


1.05 


31.7 




0.1 


4.0 


0.0 


1.0 


8 


0.5 


96 


0.0049 


0.066 


1.05 


31.7 




0.1 


5.0 


0.0 


1.0 


8 


0.5 


96 


0.0049 


0.065 


1.05 


31.7 




1.0 


1.0 


0.0 


1.0 


8 


0.5 


96 


0.0049 


0.068 


1.41 


16.0 




1.0 


2.0 


0.0 


1.0 


8 


0.5 


96 


0.0049 


0.067 


1.41 


16.0 




1.0 


4.0 


0.0 


1.0 


8 


0.5 


96 


0.0049 


0.067 


1.41 


16.0 


AD 


0.2 


0.5 


0.0 


0.5 


5 


0.5 


24 


0.0020 


0.051 


1.10 


30.8 




0.2 


0.5 


0.5 


0.5 


5 


0.5 


24 


0.0020 


0.051 


1.10 


30.8 




0.2 


0.5 


1.0 


0.5 


5 


0.5 


24 


0.0020 


0.051 


1.10 


30.8 




0.2 


0.5 


2.0 


0.5 


5 


0.5 


24 


0.0020 


0.051 


1.10 


30.8 


BD 


0.1 


1.0 


0.0 


1.0 


8 


0.5 


96 


0.0049 


0.065 


1.05 


31.7 




0.5 


1.0 


0.0 


1.0 


8 


0.5 


96 


0.0049 


0.066 


1.22 


25.6 




1.0 


1.0 


0.0 


1.0 


8 


0.5 


96 


0.0049 


0.068 


1.41 


16.0 




4.0 


1.0 


0.0 


1.0 


8 


0.5 


96 


0.0049 


0.067 


2.24 


1.9 



peak amplitude at the injection scale vi obtained from 
the spectra of velocity at the final time of simulation 
(t=12.0), the total Alfven speed Va = Vi^ox + ^oz) /Po 
and the plasma beta parameter /3 = p/pmag- Note that 
the Alfven speed Va ~ \Box\/y/po defined by the strength 
of the antiparallel component of magnetic field is the 
same for all models and equals 1.0. 

All models presented in this section were calculated 
with the grid size Ax « 0.004 corresponding to the 
resolution 256x512x256. In the process of constructing 
our numerical model we have performed number of tests 
for convergence of numerical solution during the Sweet- 
Parker stage. These studies indicated that the conver- 
gence is reached for both HLL and HLLD solvers at 
this resolution, but the HLL scheme characterizes by 
a much higher numerical dissipation. This is expected 
since the HLL solver smooths out the rotational discon- 
tinuities resulting from the presence of magnetic field. 
As will will show in § 15.11 the numerical resistivity in 
the models calculated with the HLLD scheme has the 
value of about 6.1 • 10-"*. The convergence studies signify 
that the same models calculated using the HLL scheme 
with the same resolution reveal a higher numerical re- 
sistivity of order of 8.7 • 10-^, which is over 40% larger 
than rjnum in the HLLD scheme for the same grid size 
Ace « 0.004. For models with lower resolution where the 
grid size Ax « 0.008 we obtained values of the numerical 
resistivity yy^Jf^, « 6.6 • 10"'' and t^^'^'^ « 1.12 • 10-^ for 



the HLLD and HLL schemes, respectively. Here, the dif- 
ference of numerical resistivity between both resolutions 
for the HLLD solver is insignificant, justifying the con- 
vergence of the solution, while the HLL scheme behaves 
rather poorly showing an increase of numerical dissipa- 
tion by over 75% with respect to the model with smaller 
grid size. Thus all our models presented here were calcu- 
lated using the HLLD scheme if not indicated otherwise. 

4. RECONNECTION RATE MEASURE 

We measure the reconnection rate by averaging the 
inflow velocity Vi„ divided by the Alfven speed Va over 
the inflow boundaries, i.e. 



iVn/VA) 



dxdz 



v=yii 



Va 



v=y„ 



(9) 

where S defines the area of the XZ inflow boundaries. 
Since we have two XZ boundaries, located a.t y — ymin 
and y = ymax , we need to take half of the resulting inte- 
gral. This measure works well for laminar reconnection, 
when the system is perfectly stable and where the time 
derivative of the magnetic flux is zero. In the presence 
of turbulence, however, this time derivative can fluctuate 
or the turbulence in the center of the box could affect the 
flow of the plasma. In this way we would get a flow of 
magnetic flux without the presence of reconnection. In 
order to include all effects contributing to the change of 
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Fig. 4. — Schematic 3D visualization of the reconnection rate 
evaluation. and A_ areas are defined by the sign of Bj, com- 
ponent. 

magnetic flux, we define a new more general measure of 
the reconnection rate. 

We start by considering a conserved quantity; the mag- 
netic flux First, we consider the flux contained within 
a plane inside the simulation volume (see Fig. [4]). If x is 
the direction of the reconnecting field, then we start by 
considering the time derivative of the net flux of ■ It 
is 



dt^ = - 



(j) E-dl= <j> {v X B - -qj) ■ dl 



(10) 



This equation is not exactly satisfied in our simula- 
tions unless we allow for numerical resistivity. This dis- 
crepancy can be used to derive an independent measure 
of the effective resistivity of the code, which is roughly 
consistent with value derived from the speed of laminar 
reconnection. Now we split the area of integration into 
two pieces, A+ and defined by the sign of (see 
Fig. 2]). Instead of adding two areas we subtract them, 
i.e. 



\B,\dA, 



(11) 



which we can write explicitly in terms of line integrals 
around and A_ 



dt I \B^\dA^ (p E-dl 



E-dl^ 



(12) 



: (j) sign{B^)E- dl+ J 2E-dkfac 



where k face is the line separating A+ and A- (see Fig.|4]). 
The last term describes the mutual annihilation of posi- 
tive and negative B^ along the line separating them and 
by definition, this is the reconnection rate. Note that 
this includes the motion of already reconnected flux lines 
through the plane of integration. Rather than try to 
calculate it numerically, we define the interface term as 
—2Vrec\Bx,oc\Lz, whcrc |i3s,oo| is the asymptotic abso- 
lute value of Bx, and is the width of the box. We 
can then calculate the other terms which do not involve 
trying to find the interface and the parallel component 
of the electric field. The end result, which is the new 
measure of reconnection rate, is 



1 



(j) sign{Bx) E- dl- dt J \Bx\dA 
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Fig. 5. — Evolution of the total mass M, and kinetic and mag- 



netic energies, Ekin and En 



respectively, during the Sweet- 



Parker stage. The kinetic energy i?jji„ has been amplified by a 
factor of 100 to visualize its evolution more clearly. The resistivity 
in this model is set to rju = and the guide field Sqz = 0.1. 

The electric field vx B—r]j can be further divided into 
an advection term vx B^x, a shear term vx {Byij + B^z), 
and a resistive term —rjj. With this in mind the line 
integral can be rewritten as 



sign (Bx) E- dl - 



\Bx \ {v± X x)-dl 



(14) 



- j> sign {Bx) Vx {x X Bj_) ■ dl- <j> rjj ■ dl. 



This new reconnection measure contains the time 
derivative of the absolute value of B^, and a number 
of boundary terms, such as advection of B^ across the 
boundary and the boundary integral of the resistive term 
r/j. The additional terms include all processes contribut- 
ing the time change of \Bx\. In particular, they can have 
non-zero values. 

5. RESULTS 

In this section we describe the results obtained from 
our three dimensional simulations of magnetic reconnec- 
tion in the presence of turbulence. First, we investigate 
Sweet-Parker reconnection, the stage before we inject 
turbulence. A full understanding of this stage is required 
in order to perform further analysis of reconnection in the 
presence of turbulence. 

5.1. Sweet-Parker Reconnection 



(13) 



As we described in WS.'2[ Sweet-Parker reconnection de- 
velops in our models as a result of an initial vector poten- 
tial perturbation. In order to reliably study the influence 
of turbulence on the evolution of such systems, we need 
to reach the stationary Sweet-Parker reconnection before 
we start injecting energy. 

Figure [5] shows the evolution of total mass, kinetic and 
magnetic energies until we start injecting the turbulence 
(i.e. t — 7). All quantities, after some initial adaptation, 
reach almost constant values. We remind the reader, 
that the system evolves in the presence of open boundary 
conditions, which do not guarantee perfect conservation 
of mass and total energy. Nevertheless, conservation of 
these quantities is well satisfied during the Sweet-Parker 
stage in our models. 

The reconnection rate, shown in Figure [6] (solid thick 
line representing the new reconnection rate measure), 
also confirms that we have reached a stationary state. 
Initially, the reconnection rate Vrec grows until time 
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2 J 4 b 6 

Time [tj = 1 / VJ 

Fig. 6. — Evolution of the reconnection rates, (Vin/VA) (old) and 
Vrec (new), for the same model as in Fig. [5] The reconnection rate 
grows initially until it reaches the stationary solution. Note, that 
both reconnection measures coincide during the later evolution. 

t 1.0, when it reaches a maximum value of ~ 0.05. 
Later on, it drops a bit approaching a value of 0.04. Dur- 
ing the last period of about three Alfven time units, the 
change of the reconnection rate is very small. We assume 
that these conditions guarantee a nearly steady state evo- 
lution of the system. In this plot we also show the evo- 
lution of old reconnection measure (Vm/V^) which do 
not include time variation of the flux. We see, that after 
the initial discrepancy they coincide and show the same 
value of reconnection speed. 

In order to check the reliability of our code, we estimate 
characteristic parameters of the laminar reconnection, 
such as thickness of current sheet S, and the inflow and 
outflow speeds, Vm and Vout, respectively (see Fig. [7]), 
from the profiles of the absolute value of current density 
\j\ along the Y direction and the profiles of Vx and Vy 
along the X and Y directions, respectively. We compare 
these estimates to the values obtained from the equations 
describing Sweet-Parker model. We limit this procedure 
to one model with a resolution of 256x512x256, uniform 
resistivity 77^ = 10"'^ and shared component of magnetic 
field B02 = 0.1. We analyze the snapshot obtained at 
the time t = 7 before we start introducing turbulence. 

First, using the Ohmic dissipation formula Vm — 
r]/d we estimate the total resistivity coefficient includ- 
ing the numerical contribution. We obtain a value of 
V = Vu + Vnum = Vin^ ~ 1-64 • 10^"^. It indicates that 
the numerical resistivity rjnum has a value of about 60% 
of the explicit uniform resistivity. We use this value 
of resistivity in the subsequent est imates. Next , using 
mass conservation, VmL = VoutS (|Parkeil flQSTl : I Sweet! 
|1958[) , we estimate the outflow speed Vout ■ Since the size 
of our box along X direction is L = 1.0, this gives a 
value of Vout ~ 0.394, which agrees well with the value 
K>«t = 0.395 obtained from the profile of Vj.. Since we 
already know the characteristic parameters Vout, L, and 
77, we can estimate the Lundquist number S = LVout/v 
for our models. Substituting estimated Vout, L, and r/ 
we obtain a value of S' ~ 236. As expected, our nu- 
merical models are far from the conditions observed in 
typical astrophysical objects where S ^ 10^*^ — 10^°. 
Nevertheless, we can still use them to determine their 
dependence on the characteristics of the environment. 
From Vout and S we can estimate the Sweet-Parker re- 
connection rate Vrec = VoutS~^^^ , which gives a value of 
~ 0.0251, matching almost perfectly the value of Vin fit 
from the model. These estimates indicate, that our nu- 
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Fig. 7. — Estimates of the current sheet thickness 5 (top), the 
inflow (middle) and outflow (bottom) speeds Vi„ and Vout, re- 
spectively at the time t = 7 for the model with r/u = 10^^ and 
Boz =0.1. 

merical model exhibits Sweet-Parker reconnection during 
the stationary phase. 

Before presenting our results in the presence of turbu- 
lence, we describe the topology of the flow and magnetic 
field in the laminar case. In Figure [TO] we present the 
velocity and magnetic field configuration of the steady 
state. In the left panel we show the topology of the 
velocity field as textures. The brightness of texture cor- 
responds to the amplitude of the field, while the tex- 
ture itself shows the direction of the field lines. The 
topology of the velocity field is mainly characterized by 
strong outflow regions along the midplane. This outflow 
is produced by the constant reconnection process at the 
diffusion region near the center and the ejection of the 
reconnected magnetic flux through the left and right X 
boundaries. The system is in a steady state when the 
flux, which reconnects, is counterbalanced by the incom- 
ing unreconnected flux. The inflow is much slower than 
the outflow, but its direction is still apparent in the left 
plot of Figure [TOl 

The topology of the magnetic field is shown in the mid- 
dle panel of Figure[TOl We recognize the antiparallel con- 
figuration of the magnetic lines with uniform strength 
out of the midplane region. Near the midplane, the hor- 
izontal magnetic lines are reconnected generating the Y 
component, which is ejected by a strong outflow. In addi- 
tion, we show the absolute value of current density in the 
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Fig. 8. — Evolution of total mass M, and kinetic and magnetic 
energies, Ekin and Emag, respectively for the stage with turbu- 
lence. The kinetic energy E^in is not amplified in this plot. Two 
dotted vertical lines bound the period of gradually increasing tur- 
bulence. The resistivity in this model is set to r; = 10~^ and the 
shared component of magnetic field Bqz =0.1. 

right panel of Figure [TOl We see an elongated diffusion 
region in the middle of the box, where the reconnection 
process takes place. The maximum of \j\ does not ex- 
ceed a value of 25. Evidently our boundary conditions 
allow for physically realistic inflow and outflow around a 
central current sheet. 

5.2. Effects of Turbulence 

5.2.1. Evolution in the Presence of Turbulence 

The essential part of our studies covers the effects of 
turbulence on reconnection. Our goal was to achieve a 
stationary state of Sweet-Parker reconnection, described 
in the previous subsection, and then introduce turbulence 
at a given injection scale Imj cc k~^p gradually increasing 
its strength to the desired amplitude corresponding to a 
turbulent power Pmj- We inject turbulence in a region 
surrounding the midplane and extending to the distance 
of around one quarter of the size of the box. This, natu- 
rally, limits the injection scale to kf < 3. The transition 
period, during which we increase the strength of turbu- 
lence has a length of one Alfvenic time unit in all models 
and in the presented model starts at t — 6. This means 
that from t — 7 we inject turbulence with the maximum 
power Pinj. 

In Figure [5] we present an example of the evolution of 
total mass, and kinetic and magnetic energies in a model 
with Pin = 1.0, kf = 8, and rju = 10~^. We inject 
turbulence, gradually increasing its strength from t = 6 
to t = 7. This period is marked by two dotted vertical 
lines in Figure [51 We see an increase of kinetic energy 
during this period due to the injection and after i = 7 all 
quantities saturate and do not vary much. This means 
that even in the presence of turbulence and with open 
boundary conditions our system conserves total mass and 
energy relatively well. 

In Figure [5] we show the evolution of reconnection rates 
{Vin/VA) and more advanced Vrec, both described in SQl 
In this plot we recognize an increase of both rates during 
the introduction of turbulence. After the initial period 
between t — 6 and t — 8, during which the system is 
adjusting to the new state, we see that both measures 
coincide and even though they arc fluctuating, they reach 
a stationary state characterized by faster reconnection. 
A more detailed comparison of the simple and advanced 
reconnection rate measures is presented in ^5.4\ 
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Fig. 9. — Evolution of the reconnection rates, (Vi„/VA) (old) and 
Vrec (new), for the same model as in Fig. [5] In this plot we present 
the measured rates of the Sweet-Parker reconnection Vrec,SP ^'^'^ 
during the presence of turbulence, Vrec,LV- Symbol SVrec,LV is 
the time variance. AVrec,LV is the estimated uncertainty of the 
measure. 

In Figure [9] we also show how we measure the rates 
of Sweet-Parker reconnection K-ec,SP and LV99 model 
Vrec.LV- Bccause the reconnection rates fluctuate in the 
presence of turbulence we also measure their time vari- 
ance 5Vrec.LV usiug Standard deviation. In addition to 
the time variance of Vrec, we measure their errors by 
splitting the averaging region into two subregions and 
after averaging the rates Virec and V2rec over each sub- 
region (see Fig. [5]), we take the absolute value of their 
difference AV^ec — Virec — V2rec- This difference corre- 
sponds to the error of Vrec, i.e. it is different from zero 
if the rate is not constant in time. In all further anal- 
ysis and presented plots of dependencies we use values 
estimated in this way. 

Before we discuss the main results obtained in these 
studies, we present and describe the configuration and 
topology of the velocity and magnetic fields in the pres- 
ence of turbulence. In Figure [TT] we show examples of 
XY-cuts (upper row) and XZ-cuts (lower row) through 
the box of the velocity (left panel) and magnetic field 
(middle panel) topologies with the intensities corre- 
sponding respectively to the amplitude of perpendicular 
components of velocity and magnetic field to a normal 
vector defining the plotted plane. 

The first noticeable difference compared to the Sweet- 
Parker configuration is a significant change of the velocity 
and magnetic field topologies. Velocity has a very com- 
plex and mixed structure near the midplane, since we 
constantly inject turbulence in this region (see the left 
panel in Fig. Ilip . Although the structure is very com- 
plex, most of the velocity fluctuations are perpendicular 
to the mean magnetic field. This is because we are in the 
nearly incompressible regime of turbulence (large plasma 
(3, see Tab. [1]) and most of the fluctuations propagate as 
Alfven waves along the mean magnetic field. Slow and 
fast waves, whose strengths are significantly reduced, are 
allowed to propagate in directions perpendicular to the 
mean field as well. As a result most of the turbulent ki- 
netic energy leaves the box along the magnetic lines. We 
observe, however, the efficient bending of magnetic lines 
at the midplane (see the upper middle plot in Fig. [TT|) . 
This is not result of a driving, but result of reconnection. 
In general the interface between positively and negatively 
directed magnetic lines is much more complex than in the 
case of Sweet-Parker reconnection. This complexity fa- 
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Fig. 10. — Topology and strength of the velocity field (left panel) and magnetic field (middle panel) during the Sweet-Parker reconnection 
at t = 7. The strength is calculated from the components of V and B perpendicular to the normal vector of the XY-plane. In the right 
panel we show the absolute value of current density \j\ overlapped with the magnetic vectors. The images show the XY-cut through the 
domain at Z = at time t = 7 for a model with Bqz = 0.1, rju = 10~^, rja = 0.0, and the resolution 256x512x256. 
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Fig. 11. — Topology and strength of the velocity field (left panel) and magnetic field (middle panel) in the presence of fully developed 
turbulence at time t = 12. In the right panel we show distribution of the absolute value of current density \J\ overlapped with the magnetic 
vectors. The images show the XY-cut (upper row) and XZ-cut (lower row) of the domain at the midplane of the computational box. 
Turbulence is injected with power Pi„j = 1 at scale ki„j = 8. Magnetic field reversals observed are due to magnetic reconnection rather 
than driving of turbulence, which is subAlfvenic. 
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vors creation of enhanced current density regions, where 
the local reconnection works faster since the current den- 
sity reaches higher values (see the right panel of Fig. [TT|) . 
Since we observe multiple reconnection events happening 
at the same time (compare the right panel of Fig. [Tl] to 
the Sweet-Parker case in Fig. [T0|) . the total reconnection 
rate should be significantly enhanced. 

5.2.2. Dependence on Turbulence Strength 

We ran several models (models PD in Tab. [T]) with 
varying powers of turbulence. All other parameters were 
kept the same, allowing us to estimate the dependence 
of the reconnection rate Vrec on the power of injected 
turbulence Pinj- 

Figure [12] (left panel) shows the evolution of reconnec- 
tion speed in models with turbulent power Pinj varying 
in range by more than one order of magnitude, from 0.1 
to 2.0. The evolution of Vrec reaches stationarity in a rel- 
atively short period, about one Alfven time unit, as es- 
timated from the plot. We averaged Vrec from t ~ 8.2 to 
t = 12. From this plot we see that the reconnection rates 
are relatively stable. In the right panel of Figure fT2l we 
plot the averaged reconnection speed over the strength of 
turbulence. Filled symbols represent the averaged recon- 
nection rate in the presence of turbulence. Open symbols 
represent the reconnection rate during the Sweet-Parker 
process, i.e. without turbulence. The error bars show 
the time variance of Vrec- The size of symbols indicates 
the uncertainty in our estimate of the reconnection speed 
^Vrec,LV normalized to the uncertainty in the reconnec- 
tion speed during the Sweet-Parker evolution lSVrec,SP 
using formula size — 2.0 — In AV^ec,Ly/ In AT^ec.SP- 
Meaning that if AVrec,LV is of order of AVrec,SP their 
symbols have the same size. 

Fitting a power law to the calculated points indicates 

1 /2 

that the reconnection speed Vrec scales as ~ P^^j , as 
expected from the LV99 prediction (see Eq. [4]). 

5.2.3. Dependence on Injection Scale 

We performed similar studies to determine the depen- 
dence of the reconnection speed Vrec on the scale at which 
we inject turbulence, Imj- Keeping the same power of 
turbulence for all models SD (see Tab. [T]) we inject tur- 
bulence at several scales, from kinj = 5 to kinj — 25. At 
the upper end of this range the turbulence barely broad- 
ens the Sweet-Parker current sheet. At the lower end the 
turbulent eddies are barely contained within the volume 
within which we excite turbulent motions. 

In Figure [13] (left panel) we present the results for this 
series of models. From the plot we clearly see a strong de- 
pendence of the reconnection rate on the injection scale. 
After the injection of turbulence, models with smaller 
power injection scales reach smaller values of the recon- 
nection rate Vrec- In the model with the largest injection 
scale kinj = 5, reconnection is 5 times faster than the 
Sweet-Parker rate. After averaging over time, we plot 
the dependence of Vrec on the injection scale in the right 
plot of Figure [13] 

The prediction of LV99 model is that the reconnec- 
tion rate should scale as ~ linj. In the case of a strong 
guide field we see that this relation is observed for small 
injection scales. However, when linj is comparable with 
length of the volume of driven turbulence the dependence 



becomes flatter. For a weak guide field the dependence 
is flatter from the very beginning and is better flt by 

Vrec ^ ^inj- This difference requires further investiga- 
tion. 

We see several possible sources for the discrepancy. For 
instance, the existence of a turbulent inverse cascade can 
modify the effective linj- In addition, reconnection can 
also modify the characteristics of turbulence, such as the 
power spectrum and anisotropy. We shall study this is- 
sue elsewhere, but within this work get satisfied with the 
qualitative agreement of the predictions and the numer- 
ical simulations. 

5.2.4. Dependence on Resistivity 

In the global constraint on the reconnection rate 
(Eq. [3]) derived in LV99, there is no explicit depen- 
dence on the resistivity. In order to test this, we per- 
formed another set of models (RD in Tab. [1]) in which 
we change the uniform resistivity only. We expect 
that in these models we should see at early times the 
theoretical dependence of the Sweet-Parker reconnection 

1/2 

rate {Vrec,SP ^ Vu )i ^md at late times, after we intro- 
duce turbulence, a complete lack of any dependence on 
resistivity. 

In the left panel of Figure [14] we show the evolution 
of the reconnection rate for a subset of RD models with 
Bqz — 0.1 and with r/u varying from 5 • 10""* to 5 • 10""^. 
In this plot we can recognize a strong dependence of 
the Sweet-Parker reconnection on the uniform resistiv- 
ity, while in the presence of turbulence the reconnec- 
tion speed, even though it fluctuates, reaches roughly 
the same level for all models. This indicates that, as 
expected, there is little or no dependence of the recon- 
nection rate on resistivity in the presence of turbulence. 

In the right panel of Figure [14] we show reconnection 
rates obtained from the set RD of models with a weak 
and strong guide fields, Bqz = 0.1 and Bqz = 1.0, re- 
spectively. The open symbols show the rates of laminar 
reconnection Vrec.SP- They follow precisely the theoreti- 

1/2 

cal relation Vrec.SP ~ Vu - Evolving the models further 
in the presence of turbulence allowed us to test the de- 
pendence of the reconnection rates on resistivity in the 
presence of turbulence. We plot these results using filled 
symbols. We see that there is virtually no dependence on 
r/„. Moreover, we see that this property does not depend 
on the strength of Bqz- 

These results were obtained in models with a decay- 
ing zone of resistivity near the X boundary. In order to 
exhibit the importance of properly treating the resistive 
terms at the outflow boundary conditions we show in 
Figure [15] values obtained in models with regular open 
boundary conditions without a zone of decaying resis- 
tivity. Here, the open circles show the rates of laminar 
reconnection with no turbulence. We clearly see a grow- 
ing departure of these results from the theoretical pre- 

1/2 

diction, Vrec,SP ^ Vu (solid line) for models with large 
values of uniform resistivity. For rju ^ 2 ■ 10^'^, lami- 
nar reconnection seems to be insensitive to the value of 
resistivity. This effect is purely numerical and caused 
by improper handling of boundary conditions when the 
resistive term in the induction equation starts to domi- 
nate. Open diamonds in Figure [151 show the correspond- 
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Fig. 12. — Left: Time evolution of the reconnection speed Vrec for models PD (see Table [TJ with different powers of turbulence Pinj 
(see legend). Right: The dependence of the reconnection speed Vrec on Pi„j. Error bars represent the time variance of Vrec- The size of 
symbols corresponds to the error of Vrec (the way we calculate errors is described in i|5.2ll . 




Fig. 13. — Left: Time evolution of the reconnection speed Vrec for models SD (see Table[T]l with different injection scale li„j (see legend). 
Right: The dependence of the reconnection speed Vrec on . Error bars and the size of symbols have the same meaning as in Fig. 1121 




Fig. 14. — Left: Time evolution of the reconnection speed Vrec for models RD (see Table (TJ with different resistivity coefficients r]u (see 
legend). Right: The dependence of the reconnection rate Vrec on the uniform resistivity r]u- The solid line shows the theoretical dependence 

1 /2 

in Sweet-Parker model, Vrec SP ~ Vu ■ Open symbols show Sweet-Parker reconnection rate for different models. Filled symbols show 
reconnection rates Vrec in the presence of turbulence. Dashed line shows the mean value of reconnection rate averaged over all models 
during the turbulence stage (all filled symbols). Error bars show the measure of time variance of Vrec- 
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Fig. 15. — Dependence of the reconnection rate Vrec on the 
uniform resistivity r]^ without the decaying resistivity zone. The 
solid line shows the theoretical dependence in Sweet-Parker model, 
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Vrec SP ~ Vu ■ Open circles show Sweet-Parker reconnection rate 
for different models. Open diamonds show reconnection rates Vrec 
in the presence of turbulence. Filled symbols show values corrected 
by the factor calculated from the theoretical dependence (see text). 
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Fig. 16. — Dependence of the reconnection rate Vrec on the 
anomalous resistivity rja- Two horizontal lines show the mean re- 
connection rate averaged over all plotted models during the Sweet- 
Parker stage and in the presence of turbulence (dotted and dashed 
lines, respectively). The critical value of the current density jcrit 
is set to 25.0 in all models plotted here. 

ing rates for models with turbulence. Even though ini- 
tially, for smaller values of resistivity, the reconnection 
rate seems to be independent of rju, later on, when the 
resistivity is larger, it starts to decay with ?7„. The prob- 
lem of boundaries affects the turbulent case as well. 

We can test our ability to compensate for the affect 
of large resistivity on the boundary conditions using the 
theoretical dependence Vrec.sp r^^"^ ■ If we calculate 
correction factors by taking the ratio of ?7^^^/Kec,SP, 
we can use these coefhcients to correct reconnection rate 
during the turbulent stage. The result is plotted with 
the filled symbols. The relation signifies virtually no de- 
pendence of the Yrec on r\. 

In addition to the uniform resistivity dependence, we 
have studied the dependence on anomalous effects as 
well. The results of these studies are presented in Fig- 
ure [16] and show four models with the same uniform re- 



sistivity ?7m = 5 • 10^^ and the critical current density 
jcrit = 25.0 but with the anomalous resistivity parame- 
ter r\a varying between 0.0 and 2 • 10"'^ (see models AD 
in Tab. [1]). In Figure [TBI we plot the dependence of the 
reconnection rate on the anomalous resistivity parame- 
ter r^a- We see that also the reconnection speed Vrec is 
insensitive to the value of r]a to within the variations of 
the reconnection rate in each model (see the error bars). 

In this section we shown that the results of our resis- 
tivity studies agree with the LV99 model and indicate 
no sensitivity to Ohmic resistivity in the presence of tur- 
bulence. This is important for two reasons. First, it 
supports the stochastic reconnection model proposed in 
LV99. Second, it gives reason to believe that simulations 
of reconnection and MHD turbulence in astrophysical ob- 
jects can be successfully applied to real objects, despite 
the enormous difference between the magnetic Reynolds 
numbers that can be simulated and the values that ac- 
tually obtain in stars, accretion disks, and galaxies. 

5.3. Role of the Guide Field 

The robustness of the LV99 model is determined by 
its application to any configuration of the reconnecting 
field. Usually textbooks describe the Sweet-Parker model 
as having opossitely directed field lines, which undergo 
the reconnection process in the diffusion region and are 
perfectly antiparallel, since this model is strictly two di- 
mensional. This situation is very particular. In reality, 
even in the absence of turbulence, magnetic field lines can 
enter the diffusion region at different angles a (see right 
panel of Fig.[3|), or in other words, an uniform component 
parallel to the reconnecting field could be present. This 
component is called a guide field and in our simulations 
is determined by the value of Bqz- 

Direct simulations of the LV99 model with a strong 
guide field of order of the reconnecting component 
strength encounter problems when the Z boundaries are 
open. In our picture, the simulation domain is a box 
embedded in a large scale configuration of magnetic field 
around the interface of two volumes characterized by the 
different direction of magnetic lines. In this picture, mag- 
netic lines crossing our domain extend to infinity. This 
means that any force acting on a field line in the domain 
will feel the tension of that line, so the general configura- 
tion of the field line will not change. It only encounters 
a small local perturbation. When the open boundary 
conditions are applied, the force acting on a magnetic 
line is not bounded by the magnetic tension anymore, 
and in the presence of shear flows generated by turbu- 
lence, magnetic line ends attached to the boundaries can 
slip along them changing the global topology of mag- 
netic field. Moreover, the change of direction of the field 
transforms Bz into component, which suppresses the 
inflow of fresh unreconnected flux. This is an artificial 
and undesirable situation. 

There are two solutions to this problem. First, we 
can construct a special type of open boundary condi- 
tions, which fixes magnetic lines at some distance far 
from the box. Then, at the boundary, we would have to 
extrapolate the magnetic field lines and calculate the re- 
quired tension from this condition. This approach, how- 
ever, seems to be complex and requires additional condi- 
tions assumed a priori at the boundaries. Thus, we have 
chosen another approach by applying periodic boundary 
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Fig. 17. — Demonstration of the importance of boundary con- 
dition along Z-direction. Upper: Comparison of two models with 
the same conditions but different types of Z boundaries for a weak 
Bqz- Lower: Comparison of two models with the same conditions 
but different types of Z boundaries for a strong Bq^ . 

conditions along the Z direction. This is much simpler 
and does not limit our model, i.e. it does not influence 
the inflow and outflow of the magnetic flux. 

In FigurefTTlwe show the comparison of Vrec for models 
with a weak (top) and strong (bottom) guide fields Bqz. 
Each plot contains two otherwise identical models, but 
one with open boundaries along Z direction and another 
with periodic boundaries. As we see, in the case of weak 
Bqz there is no difference if we apply open or periodic 
boundaries. In both models, mean values of reconnec- 
tion speed during the presence of turbulence are almost 
identical (see black and grey dashed lines in top panel of 
Fig. fT7|) . However, once we start increasing the strength 
of the guide field, our choice of the boundary conditions 
starts to be important. In the case of open boundaries 
along the Z direction, the mean reconnection rate is sup- 
pressed by the strong turbulent shear of B^, which is 
clearly seen in the bottom plot of Figure [T71 Applying 
periodic boundaries along the Z direction we can restore 
the fast reconnection rate, but the system needs more 
time to adjust to a steady state. This is caused by the 
fact that after starting injection of turbulence, part of the 
injected kinetic energy is converted into magnetic energy, 
thus the relaxation of the system takes place more slowly 
than in the case of weak B^. This demonstration shows 
that the choice of applied boundary conditions along the 
Z direction is important and has a crucial infiuence on 
the reconnection rate, at least in the case of strong mag- 
netic guide field. 

In Figure [18] we show four models (set BD in Tab. [ij 
with different strengths of the guide field Bq^. This plot 
shows that the reconnection speed is comparable over 
a broad range of Bqz (and the angle a), which means 
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Fig. 18. — Dependence of the reconnection rate Vrec on the 
strength of the guide field Bqz • 

that reconnection will work with similar efficiency for 
any configuration of magnetic field in which the lines of 
the magnetic field enter the diffusion region at different 
angles. The studied models with 0.1 < Bq^ < 1.0 cover 
range of a from 0° to 90°. For a given power input a 
larger guide field is expected to correspond to a smaller 
amount of field wandering and therefore a smaller recon- 
nection speed (see the model with Bqz = 4.0 correspond- 
ing to a ~ 152°). We, however, do not see a prominent 
dependence on the shared component of magnetic field, 
which calls for further studies of the effects of the strong 
guide field. 

5.4. Comparison of the Old and New Measures 

In Figure [19] we show a comparison of the reconnec- 
tion rate obtained using the old and new measures as a 
function of the power of turbulence (left plot) and the 
injection scale (right plot). In the case of the turbulent 
power dependence we used two different sets of models 
for each method. For the old method of Vrec estima- 
tion the guide field Bz = 0.2, and the uniform resistivity 
rj = 5 ■ 10~^. The guide field in new models is set to 0.1. 
The injection scales are kinj = 5 and 8, for the old and 
new models, respectively. Dashed line shows the LV99 
dependence Vrec ~ V^. 

The old measure models were calculated with a more 
dissipative scheme based on the HLL solver instead of 
HLLD, with a higher value of numerical resistivity. Dur- 
ing the Sweet-Parker stage we note some difference in the 
values of reconnection rate between two sets of models 
resulting from the different values of uniform resistivity 
rju in the models and other parameters, and different type 
of schemes used to evolve the numerical models. After 
the injection of turbulence, as we shown in Figures [Tl] 
and [TH] during the stage of the presence of turbulence, 
the reconnection rate is insensitive to the value of resis- 
tivity. In Figure [T9I we see that both relations, using the 
old and new methods, show the same dependence, even 
though they were fitted to two different sets of models. 
This means that the reconnection rate dependence on 
the injected power is not sensitive to the strength of the 
guide field Bz, the injection scale or the value of uniform 
resistivity and the numerical diffusivity of the method 
used in solving the MHD equations. 

In the right plot of Figure[T9]we see a comparison of the 
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Fig. 19. — Comparison of the old and new reconnection rate estimates. In the left plot we show the dependence of the reconnectlon 
rate Vrec on the power of turbulence Pi„j (diamonds) using both methods. The corresponding Sweet-Parker rates, without the presence 



of turbulence, are shown using x-symbols. In the right plot wc show the dependence of the reconnection rate Vr, 
linj (diamonds). Again, the Sweet-Parker rates are shown using x-symbols. 

old and new measures of Vrec as a function of the injec- 
tion scale linj- Similar to the dependence on the power 
of turbulence, this dependence does not change and in 

both cases is Vrec ~ This confirms again that the 

dependence of the reconnection rate Vrec on the proper- 
ties of the turbulence do not change its character. The 
only difference seen in both relations is the amplitude of 
the reconnection rate. This difference comes from the 
fact, that in the old measure models we used weaker tur- 



on the injection scale 



bulence {Pinj — 0.5 vs. Pi, 



1.0). Another difference 



could result from the fact that we include more terms in 
the elaboration of the new reconnection rate, such as the 
boundary shear, the boundary resistive terms, the time 
derivative of the absolute value of iBa^j. These terms can 
have non-zero time averages, although we expect them 
to have at least stationary values. 

The results showed here confirms that the new and old 
methods of estimating the reconnection rate reveal the 
same dependence on the power of turbulence and injec- 
tion scale. It signifies that the old measure {Vin/VA) and 
new one Vrec are both adequate in determining the speed 
of reconnection and more importantly, that they measure 
precisely the reconnection rate in our system excluding 
all other processes, which could contribute to the plasma 
flow, but not to the reconnection process itself. It is 
also important that these results are independent of the 
numerical scheme or the initial conditions. 

5.5. N on- Reconnection Case 

The results presented so far indicate an agreement with 
the LV99 model. We see that the presence of turbulence 
near current sheets enhances reconnection significantly. 
However, since our testing is limited to numerical mod- 
els, which have explicitly defined boundaries, our results 
may be sensitive to interactions between turbulence and 
the boundary conditions implemented in the numerical 
model. These interactions could stimulate a flow through 
the box resulting in a false, non-zero reconnection rate 
without the presence of reconnection itself. Moreover, 
the measures of reconnection rate that we applied here 
may indicate the existence of other processes contribut- 
ing to the plasma flow and related to the turbulence itself 
rather than reconnection. In order to test the reliability 
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Fig. 20. — Comparison of the reconnection rate Vrec for two 



models with the same set of parameters, Pi„ 



1.0, k 



f 



Bo. 



0.1, rju = 10~ , but different initial configurations of magnetic 
field: reversed (grey) and uniform Bx (black). Both models 
were run with the same resolution 128x256x128. 

of our numerical model and the methods of reconnection 
rate estimation we present below a comparison of two 
models with the same set of parameters, but with differ- 
ent initial configurations of magnetic field, namely, one 
with the antiparallel as used in all numerical models 
presented here and another with the uniform where is 
no laminar reconnection before we start injecting turbu- 
lence. From the latter model one would obviously expect 
the lack of reconnection and zero reconnection rate Vrec ■ 

In Figure [20] we show the evolution of reconnection 
rate Vrec for two models with the same set of parame- 
ters, Pinj = 1.0, kinj = 8, Boz = 0.1, and ?7„ = 10"^. 
Models were started from different initial configurations 
of magnetic field. The initial configuration of B in the 
first model is taken from the steady-state laminar recon- 
nection of Sweet-Parker type. In the second model we 
set uniform magnetic field B^ = 1.0 initially. Starting 
the evolution from these conditions we inject turbulence 
in the same way in both models, increasing its strength 
gradually through the period of one Alfvenic time unit. 
Then, turbulence is injected with the power of constant 
value of Pint = 1-0. 

Figure [20] shows a clear difference between the recon- 
nection rates in both models, i.e., in the case of initial 
Sweet-Parker configuration further evolution of Vrec is 
similar to the results presented in the previous sections. 
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The reconnection rate grows from a value corresponding 
to the Sweet-Parker rate of about 0.02 to higher val- 
ues when it saturates at time t « 1.4 at the level of 
Vrec ~ 0.07. Comparing this to the evolution of Vrec in 
the model with the initial uniform field we see only small 
fluctuations of reconnection rate around zero, what indi- 
cates the absence of global reconnection in the system. 

Results presented in this subsection indicate essentially 
a lack of reconnection in the case of uniform initial field, 
what means that the enhancement of reconnection rate 
in the LV99 model is resulting from an infiucncc of tur- 
bulence on the laminar reconnection process and not a 
result of turbulence itself. 



6. DISCUSSION 

6.1. Relation to Other Studies of Turbulent 
Reconnection 

Our simulations presented in this paper test the LV99 
model of reconnection. At the same time, the notion that 
magnetic field stochasticity might affect current sheet 
structures is not unprecedented. For instance, some pa- 
pers have concentrated on the effects that turbul ence in- 
duces on the microphysical level. In particular, iSpeiseH 
([1970) showed that in coUisionless plasmas the electron 
collision time should be replaced with the electron reten- 
tion t ime in the current sheet. Also iJacobson fc Mosed 
(| 19841 ) proposed that the current diffusivity should be 
modified to include the diffusion of electrons across the 
mean field due to small scale stochasticity. Our simula- 
tions are performed at the MHD level and do not take 
these effects into account^. 

" Hyper-resistivity]] (jStraussI 

Bhattacharice & Hameiril 

Hameiri fc Bhattacharied Il987t [Diamond fc Malkovl 



1986 
1986 



2001 " is a more subtle attempt to derive fast reconnec- 



tion from turbulence within the context of mean-field 
resistive MHD. The form of the parallel electric field 
can be derived from magnetic helicity conservation. 
Integrating by parts one obtains a term which looks 
like an effective resistivity proportional to the magnetic 
helicity current. There are several assumptions implicit 
in this derivation. The most important objection to this 
approach is that by adopting a mean-field approxima- 
tion, one is already assuming some sort of small-scale 
smearing effect, equivalent to fast reconnection. Fur- 
thermore, the integration by parts involves assuming a 
large scale magnetic helicity flux through the boundaries 
of the exact form required to drive fast reconnection. As 
we are not aware of any attempts to test the idea of 
"hyper-resistivity" numerically, we shall not discuss this 
c oncept further. 

IStraussI (|1988( ) examined the enhancement of recon- 
nection through the effect of tearing mode instabilities 
within current sheets. However, the resulting reconnec- 
tion speed enhancement is roughly what one would ex- 
pect based simply on the broadening of the current sh eets 
due to internal mixing. IShibata fc Tanumal (|2001h ex- 
tended the concept suggesting that tearing may result 

* These effects will usually be small compared to effect of a broad 
outflow zone containing both plasma and ejected shared magnetic 
flux. Moreover, while both of these effects will affect reconnection 
rates, they are not sufficient to produce reconnection speeds com- 
parable to the Alfven speed in most astrophysical environments. 



in fractal re c onnec tion taking place on very small scales. 
IWaelbroeci] (|1989r ) considered not the tearing mode, but 
the resistive kink mode to accelerate reconnection. 

We do not expect either tearing or kink modes 
to allow us to evade the constraints on the global 
plasma flow that leads to slow reconnection speeds^. 
The effects of tearing modes may be important for 
the initiation of turbulence in the situ ation when the 
initia l configuration is lam i nar ( se e [Dahlburg et al.l 
1991 [Pahlburg fc Karoenl [1991 iDahlburd 
Ferraro fc Rogers! 12004 IL azarian fc Vishniad l2008f ). 



We plan to study this effect elsewhere. However, once 
turbulence is initiated reconnection should proceed 
independently of the tearing mode. Straus' idea is 
closely related to recent attempts to explain the exis- 
tence of the observed thick reconnection region in the 
conte xt of numerical models of c oUisionless reconnec- 
tion (jCiaravella fc Ravmondl l2008f) and p redictions of 
partic le acceleration within these regions (jDrake et al.l 
l2f)06al) . Here we do not appeal to coUisionless effects, 
but demonstrate both thick reconnection regions and 
fast reconnection in the presence of turbulence. 

The close st study to our s was done 
by J vlattha eus fc LamkinI (|1985f) (see also 
iMatthaeus fc LamkinI 119861 ). The authors studied 
2D magnetic reconnection in the presence of external 
turbulence. An enhancement of the reconnection 
rate was reported, but the numerical setup precluded 
the calculation of a long term average reconnection 
rate. A more recen t study along the approach in 
IMatthaeus fc LamkinI |l985) is one in IWatson et al.l 
( 20071) . where the effects of small scale turbulence on 2D 
reconnection were studied and no significant effects of 
turbulence on reconnection were reported for the setup 
chosen by the authors. The major differences from the 
present study stem from the fact that we test a 3D 
model of reconnection, as the LV99 depends on effects, 
e.g. field wandering, that happen only in 3D. Thus 
within the present study we do not address the question 
of what is happening to 2D reconnection in the presence 
of turbulence, appealing instead to the 3D nature of 
astrophysical reality (see more in i|6.4p . 

Finally, we would like to stress that our calculations do 
not attempt to validate the concept of turbulent diffusiv- 
ity frequently emp loyed in t he mean field dynamo (see 
iMoffat 1978:' Black man fc Fi eld 2008). This concept ap- 
peals to the alleged ability of hydrodynamic motions to 
mix magnetic lines on very small scales, as if the magnetic 
field were a passive scalar. Naturally, this assumption is 
wrong for any astrophysical field of dynamical impor- 
tance. In our calculations, on the contrary, the magnetic 
field energy exceeds the kinetic energy and thus, the tur- 
bulence is weak. Instead of mixing magnetic lines, as 
in the case of "magnetic turbulent diffusivity" creating 
reversals of magnetic fields on the smallest scales, our 

^ Sometimes the growth rate of the resistive modes is erroneously 
identified with the reconnection speeds. This does not account for 
the global outfiow of the matter constraint that should be satisfied 
for reconnection of magnetic field. When this is done, the results 
are different. For instance, LV99 equated the shearing rates arising 
from the global outfiow and the rate of the tearing mode growth at 
the scale of the reconnection region. They obtained the reconnec- 
tion speed VaR^^^ t which was faster than the Sweet-Parker rate, 
but incredibly slow for most astrophysical systems. 
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driving is only able to induce wandering of the magnetic 
field lines within a well-defined mean direction of the 
magnetic field. 

In general, it is important not to confuse stochastic 
reconnection with the ill motivate d concept of turbu- 
lent difFusivity. For instance, Kim fc Diamond ()2001l ) 
addressed the problem of stochastic reconnection by cal- 
culating the turbulent diffusion rate for magnetic flux 
inside a current sheet. They obtained similar turbu- 
lent diffusion rates for both two dimensional and reduced 
three dimensional MHD (2D with the Z-components of 
velocity and magnetic field). In both cases the presence 
of turbulence had a negligible effect on the fiux trans- 
port. The authors pointed out that this would prevent 
the anomalous transport of magnetic flux within the cur- 
rent sheet and concluded that both 2D and 3D stochas- 
tic reconnection proceed at the Sweet-Parker rate, even 
if the individual small scale reconnection events happen 
quickly, which apparently contradicts our calculations in 
the paper. 

We would like to stress that the turbulent diffusion 
rates within the current sheet are irrelevant to the pro- 
cess of stochastic reconnection (see discussion in LVC04) . 
The basic claim in LV99 is that the realistic magnetic 
field topologies allow multiple connections between the 
current sheet and the exterior environment, which would 
persist even if the stochastic magnetic field lines are sta- 
tionary ("frozen in time") before reconnection happens. 
This leads to global outflow constraints which are weak 
and do not depend on the properties of the current sheet. 
In particular, the analysis in LV99 assumed that the cur- 
rent sheet thickness is determined purely by Ohmic dissi- 
pation and that turbulent diffusion of the magnetic field 
is negligible inside and outside of the current sheet. 

6.2. 2D versus 3D reconnection 

The fact that our study is in 3D is essential, as the 
LV99 model is intrinsically three dimensional. The gen- 
eral picture is of tangled field lines, with reconnection 
taking place via a series of " Y-points" or modified Sweet- 
Parker sheets distributed in some fractal way throughout 
the turbulence. A large scale Sweet-Parker sheet will be 
replaced by a more fractured surface, but the current 
sheets will occupy a vanishingly small fraction of the to- 
tal volume and the field reversal will remain relatively 
well localized. The model predicts that the reconnection 
speed would be approximately equal to the strong tur- 
bulent velocity with a modest dependence on the ratio 
of the eddy length to the current sheet length. There 
should be no dependence on resistivity. The major re- 
sults contained in our figures showing the dependence of 
the reconnection speed on resistivity, input power and 
input scale agree with the quantitative predictions of the 
LV99 model. We are not aware of any competing models 
to compare our simulations with. 

In the absence of quantitative model to be tested, simu- 
lations aimed at studying reconnection speed have been 
done in 2D. This allowed to achieve higher resolutions 
(compared to that contemporary available in 3D) but 
substantially constrained magn etic field dynamics. For 
instance, the 2D study in Mat thaeus fc LamkinI (jl986l ) 
stresses the importance of turbulence for modifying the 
character of magnetic reconnection and specifics heat- 
ing and transport as the effect of particular significance. 



as well as formation of Petscheck-type "X-points" in 2D 
turbulence. 

As we mentioned earlier, LV99 appeals to "Y-points" 
and heating do es play a role in it. Our simulations are 
also isothermal. iKim fc Diamondl (|2001l ) showed that the 
transport of magnetic flux is not enhanced to the recon- 
nection zone. At the same time, field wandering is the es- 
sential feature of 3D reconnection, w hich is absent in the 
2D c ase. In our forthcoming paper (jKulpa-Dvbel et ahl 
I2OO9I) we shall demonstrate the difference of 2D and 3D 
reconnection using direct numerical simulations^*^. 

6.3. Reconnection in collisionless and collisional plasma 

The LV99 model of reconnection is applicable to the 
collisional medium. For instance, it is applicable to the 
interstellar medium (known to be both turbulent and 
magnetized), for which one cannot apply the Hall-MHD 
reconnection (lYamadall2007h . This is a relieving news for 
the int erstellar medium, star formation and Solar s imula- 
tions (jOstriker et al.ll2001t |McKee fc Ostriker 2007), as 
numerical reconnection in MHD codes is fast. Note, that 
the requirement of being collisionless in terms of mag- 
netic reconnection is different from the usual use of the 
term in astrophysics. For instance, for the diffuse ISM 
the collisionallity parameter is lOcT,,, where ujc is cyclotron 
frequency of electrons and is the collisional time for an 
electron. However, for Hall-MHD reconnection to be ap- 
plicable the criterion is different. It is required that the 
Sweet-Parker current sheet 6sp width be less than the 
ion incrtial length di. Thus the "reconnection criterion 
for media to be coUisionless" is (L/di)^^^ /{uJcTe) < 1, 
which presents a much severe constraint on the possible 
rate of collisions. As a result, magnetic reconnection hap- 
pens mediated by the Hall-MHD only if the extend of the 
contact region L (see Fig. [T]) does not exceed 10^^ cm. 
Magnetic fields in the ISM should interact over much 
larger scales. 

Astrophysical environments also contain media to 
which the Hall-MHD reconnection is applicable, e.g. So- 
lar corona, interplanetary medium. Is LV99 model ap- 
plicable to such environments? The answer is positive if 
the level of turbulence is high enough. Indeed, the re- 
connection on microscales can happen fast, i.e. in the 
Hall-MHD fashion. However, this may not change the 
global reconnection rate. Indeed, the LV99 model shows 
that even with relatively slow Sweet-Parker reconnection 
at microscales the global reconnection is limited not by 
Ohmic resistivity, but the rate of magnetic field wonder- 
ing. In fact, anomalous resistivity that we use in some 
of our numerical runs is a proxy for the Hall-MHD effect 
on reconnection and we do not see any dependence of 
the reconnection rate on the anomalous resistivity (see 
Fig. [TH). We believe that the Hall-MHD local recon- 
nection of magnetic fields is taking place interplanetary 
medium, which is being tested by local in-situ measure- 
ments, while the global reconnection rates are determined 
by magnetic field wandering as prescribed in LV99. Fu- 
ture dedicated numerical experiments and space mea- 
surements should test this idea. 

6.4. Reconnection in Partially Ionized Gas 

To our best knowledge no quantitative predictions are cur- 
rently available for the 2D reconnection. 
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The LV99 model of reconnection deals with either fully 
ionized plasma or with plasmas where the neutrals con- 
stitute less than 10% of species. For higher percentage 
of the neutrals, the MHD-type turbulent cascade is trun- 
cated by the viscosity of neutrals. How does reconnection 
happen in a partially ionized gas? 

A partially ionized plasma fills a substantial volume 
within our galaxy and the earlier stages of star forma- 
tion take place in a largely neutral medium. This mo- 
tivates our study of the effect of neutrals on reconnec- 
tion. The role of ion-neutral collisions is not trivial. On 
one hand, they may truncate the turbulent cascade, re- 
ducing the small scale stochasticity and decreasing the 
reconnection speed. On the other hand, the ability of 
neutrals to diffuse perpendicular to the magnetic lines 
allows for a broader particle outflow and enhances re- 
connection rates. 

Reconnection in partially ionized gases before the in- 
troduction of t he LV99 model looked hope lessly slow. 
For instance, in i Vishniac fc LazarianI ()1999l . henceforth 
VL99) we studied the diffusion of neutrals far from 
the reconnection zone assuming the anti-parallel mag- 
netic lines. The ambipolar reconnection rates obtained 
in VL99, although large compared to the Sweet-Parker 
model, are insufficient either for fast dynamo models or 
for the ejection of magnetic flux prior to star formation. 
In fact, the increase in the reconnection speed stemmed 
entirely from the compression of ions in the current sheet, 
with the consequent enhancement of both recombination 
and ohmic dissipation. This effect is small unless the re- 
connecting magnetic fie ld lines are almost exact ly anti- 
parallel (VL99, see also lHeitsch fc Zweibef2003allbi ').Anv 
dynamically significant shared field component will pre- 
vent noticeable plasma compression in the current sheet, 
and lead to speeds practically indistinguishable from the 
standard Sweet-Parker result. Since generic reconnection 
regions will have a shared field component of the same 
order as the reversing component, the implication is that 
reconnection and ambipolar diffusion do not change re- 
connection speeds significantly. 

LCV04 presented a model of turbulence in a partially 
ionized gas. This model agrees well with numerical sim- 
ulations available as the limiting case which can be char- 
acterized by one fluid with a high Prandtl number (see 
ICho fc LazariMl[200l [200l . Using this model LVC04 
described field wandering, which is the core of the LV99 
model of reconnection. They showed that the magnetic 
reconnection proceeds fast, both in the diffuse interstel- 
lar and molecular cloud partially ionized gas. 

In the paper above we employed one-fluid code and 
therefore we could not test the predictions in LVC04. 
However, LVC04 model at its core is based on the ideas 
presented in LV99. Therefore, our successful testing of 
LV99 model provides us with optimism in relation to the 
LVC04 model. The latter should be tested with a two 
fluid code capitalizing on our present experience with 
the boundary conditions, turbulent driving and the ways 
of measuring of reconnection rate obtained in this pa- 
per. Needless to say, the numerical confirmation of fast 
reconnection in partially ionized gas would open wide 
prospects for quantifying processes taking place in molec- 
ular clouds, including the eni gmatic processes of star for- 
mation. For instance. IShu et al.. (.200 6) showed that mag- 
netic field is removed from star-forming cores faster than 



allowed by the standard ambi pola r diffusion s cenari o (see 
iTassis fc 'Mouschoviasl[2005al lHl. IShu et al.l ()2007f ) pro- 
posed a mechanism that required efficient reconnection 
of magnetic loops based on the "hyper-resistivity" con- 
cept. As we discussed earlier, this concept is not self- 
consistent and is therefore problematic. Fast magnetic 
reconnection in the partially ionized gas could serve a 
similar purpose. 

6.5. The Current Simulations and Remaining Problems 

In this paper we have performed extensive numerical 
testing of the LV99 model of magnetic reconnection. A 
advantage of this model is that it is very generic and 
does not appeal to any particular plasma effects and/or 
initial boundary conditions. In many respects, it is a nat- 
ural generalization of the Sweet-Parker model to the case 
when turbulence is present. It does not require the cur- 
rent sheet to open up, as is required in the Petscheck re- 
connection model. The only requirement is the presence 
of field line wandering, which is a well documented aspect 
of magnetic turbulence (see LVC04 and ref. therein). 

Our results indicate that the reconnection speed in- 
creases in the presence of weak turbulence in a manner 
consistent with the LV99 predictions. In particular, the 
reconnection speed Vrec (Eq. [3]) shows a strong depen- 
dence on the characteristics of the turbulence, such as 
the rate of energy dissipation and the scale of injection. 
It is very important that we observe no explicit depen- 
dence on the fluid resistivity, which is the requirement 
for fast magnetic reconnection in astrophysical environ- 
ments. 

To make progress in simulating the weakly stochastic 
reconnection, we had to develop a new numerical set up, 
which included both driving turbulence over a part of the 
box, inflow conditions for the incoming flux and outflow 
conditions for the reconnected flux. This enabled us to 
study both the transitional regimes of reconnection and 
the steady state one. In comparison, simulations in the 
periodic box would not be instructive, as stochastic re- 
connection is dynamic and does not reach a stationary 
state in one crossing time. 

In addition, we have developed a new way of measur- 
ing the reconnected flux (see §|4]). This was required as 
worries were associated with the accuracy of the "intu- 
itive measure" , which is the rate with which the unrecon- 
nected flux enters the computational box at the inflow 
boundaries. Potentially, the "intuitive measure" could 
overestimate the reconnected flux if the unreconnected 
magnetic lines escape through the leaky outflow bound- 
aries of the computation box. The new measure accounts 
for such a loss of magnetic flux. 

To our satisfaction, we found that the new and old 
measures provide consistent results both in the case of 
no turbulence and with turbulence, which means that the 
effect of the outflow of unreconnected field was subdom- 
inant for the situations when the magnetic fluxes were 
intersecting at an appreciable angle. For the limiting 
case, however, when the large scale reconnection was ab- 
sent, i.e. when the corresponding angle between the two 
fluxes was zero and without turbulence magnetic fields 
were all parallel, the "intuitive measure" was a non-zero 
value, while the new measure of reconnection gave zero 
reconnection rate. This shows the new measure is more 
reliable. 
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Our results confirm that the reconnection of magnetic 
field in the presence of turbulence is independent of 
Ohmic resistivity, allowing for fast reconnection. In ad- 
dition, we showed that the reconnection rate does not 
depend on the anomalous resistivity, which proves that 
the LV99 model of reconnection is also applicable to coUi- 
sionless plasmas and may dominate reconnection in such 
systems if the lev el of turbulence is efficiently high. Note, 
that according to lLazarian fc Vishniad (|2008f ) turbulence 
in a system may be successfully generated by the re- 
lease of magnetic energy through reconnection. Thus, we 
anticipate that the initial reconnection in a coUisionless 
plasma may be induced through the Hall reconnection, 
but later the resulting MHD turbulence will determine 
the reconnection rate. 

Apart from the prediction in LV99 that the reconnec- 
tion rate should be independent of resistivity, we have 
also tested the quantitative scalings of the reconnection 
rates presented in LV99. In particular, we tested the de- 
pendence of the reconnection rate on the turbulent injec- 
tion power. Our numerical results confirm the LV99 pre- 
diction Vrec ^ Vj^ which, as we explained earlier, trans- 
lates into Vrec Plnj (sce Eq. 13]). 

We do not expect to see the dependencies of the re- 
connection rate on viscosity until the Prandtl number is 
less than unity since the cascade can proceed to scales 
smaller than the current sheet thickness affecting the rate 
of local reconnection. This is the case we are studying 
here. For large Prandtl numbers the reconnection may 
be suppressed as magnetic field becomes more laminar at 
sufficiently large scales. However, two main factors of the 
LV99 model contributing to the enhancement of the re- 
connection rate in the presence of turbulence are that the 
turbulence allows for multiple reconnection events act- 
ing simultaneously and secondly, the turbulence should 
broaden the outflow region to the widths correspond- 
ing to the scales much larger than the dissipation scale 
and comparable to the injection scale allowing to remove 
more reconnected flux. While the flrst factor could be af- 
fected by the large viscosity, the second is still valid if the 
injection scale is much larger than the dissipation scale. 
The case when the effective Prandtl numbers are larger 
than unity due to the presence of neutrals is discussed in 
iLazarian ct al. ( 2003)- We plan to study the dependence 
on the Prandtl number in our next publication. 

The LV99 model predicts that the reconnection speed 
should grow with the injection scale as Vrec ~ knj ■ How- 
ever, some of our numerical tests show a weaker depen- 
dence, Vrec ^ liij- The origin of this discrepancy is 
unclear. Our model of how weak isotropic turbulence 
evolves cascades into strong turbulence may be too naive. 
In particular, the approximation that the parallel wave- 
length is unchanged during this process may be sim- 
plistic. In addition, we may simply not have sufflcient 
dynamic range between the dimensions of the subvol- 
ume where turbulence is excited and the thickness of the 
Sweet-Parker layer. Corrections at either end could re- 
sult in the appearance of a shallower dependence on the 
driving scale. The true nature of this discrepancy should 
be explored in future work. 

Even though these numerical simulations allow us to 
study reconnection in the presence of turbulence for a 
limited range of magnetic Reynolds numbers (in this pa- 



per Rjn < 10'^), the results provide good testing of the 
relations derived by LV99. The strong dependence of Vrec 
on the injection scale, when scaled to the real conditions 
of the interstellar medium, shows a dramatic enhance- 
ment of the reconnection speed, even in the situation of 
an almost perfectly frozen magnetic field in the medium. 
This allows for fast reconnection with the characteristic 
time comparable to one Alfven time unit. 

6.6. Implications of the LV99 Model 

Reconnection is one of the most fundamental processes 
involving magnetic fields in conducting fiuids or plasmas. 
Therefore, the identification of a robust process responsi- 
ble for reconnection has many astrophysically important 
consequences. Below we shall mention a few selected im- 
plications o f the successful valid a tion o f the LV99 model. 

Following IZweibel fc Yamadal (|2009D we note that so- 
lar flares in spired much of the earlier research on recon- 
necti on (see lPneuma 3 ll981t iBastian et a!] ' 1998': H udsonl 
|2008( ). As the plasma involved is substantially rar- 
efied, the restrictive conditions for the coUisionless re- 
connection are sati sfied in this particular environment. 
ICassak et al.l ()2005| ) stated that bistable Hall reconnec- 
tion can be important in this case. Stochastic reconnec- 
tion, as we discuss below, provides an alternative expla- 
nation. 

Indeed, an important prediction of the LV99 model is 
related to the reconnection instability that arises in the 
situation when the initial structure of the flux prior to re- 
connection is laminar. Reconnection at the Sweet-Parker 
rate is negligible. This allows magnetic flux to accumu- 
late. However, when the degree of stochasticity exceeds 
a threshold value, the reconnection itself should excite 
more turbule nce, creating a positive feedback resulting 
in a fiare (see ILazarian fc V ishniac 2008) . The instabil- 
ity is a generic property of laminar field reconnection in 
both coUisionless and coUisional environments. Referring 
to the Sun, one may speculate that the difference between 
gradual and eruptive flares arises from the original state 
of magnetic fleld prior to the fiare. In the case when the 
magnetic field is sufficiently turbulent the accumulation 
of magnetic flux does not happen and the flare is gradual. 
Similarly, the observed spatial spread of energy release 
during solar flares may be due to the spread of the re- 
gion of turbulent fields once reconnection is initiated at 
one place. Further research is necessary for testing these 
ideas. 

In the Sweet-Parker model reconnection can accelerate 
charged particles, e.g. due to the electric field in the re- 
connection region (see iLitvT nenko 2003). However, the 
speed of Sweet-Parker reconnection is negligible for most 
astrophysical environments, thus the transfer of energy 
from the magnetic field to particles is absolutely negligi- 
ble, if reconnection follows the Sweet-Parker predictions. 

It is interesting to notice that the first-order Fermi 
acceleration process is intrinsic to the LV99 model of 
reconnection. Consider a particle entrained on a re- 
connected magnetic field line (see FigH]). This parti- 
cle may bounce back and forth between magnetic mir- 
rors formed by oppositely directed magnetic fiuxes mov- 
ing towards each other with the velocity Vr. Each 
bounce will increase the energy of a particle in a way 
consistent with the requirements of the first-order Fermi 
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proce ss " (Ide Gouveia dal Pino fc LazarianI l2001l |2003| . 
[2005tlLazariaijl200fil ^. This is in contrast to the second- 
order Fermi acceleration that is frequently discussed in 
terms of acceler a.ting particles by tu rbulence generated 
by reconnection ()La Rosa et al.l[2006l ). 

An interesting property of this acceleration mech- 
anism is that it is potentially testable observation- 
ally, since the resulting spectrum of accelerated par- 
ticles is different from that arising from a shock. 
Ide Gouveia dal Pino fc LazarianI (|2001l . [2005) used this 
mechanism of particle acceleration to explain the syn- 
chrotron power-law spectrum arising from the flares of 
the microquasar GRS 1915-1-105. The mechanism is sim- 
ilar to t he accelerat i on mec h anism that was disc ussed 
later by iDrake et"aI1 (|2006bD . iDrake et all (|2006bD con- 
sidered the acceleration of electrons and, similarly, to the 
iMatthaeus et al.l ()1984D . assumed that the acceleration 
happens within 2D contracting loops. 

To enable sustainable dynamo action and, for exam- 
ple, generate a galactic magnetic field, it is necessary to 
reconnect and rearrange magnetic flux on a scale simi- 
lar to a galactic disc thickness within roughly a galactic 
eddy turnover time (~ 10® years). This implies that re- 
connection must occur at a substantial fraction of the 
Alfven velocity. The preceding arguments indicate that 
such reconnection velocities should be attainable if we 
allow for a realistic magnetic field structur e, one that in- 
clude s both random and regular fields (see iHanasz et al.l 
|2004 e.g.). This does solve one part of the problem of 
dynamo. The ot her part is related to magnetic helicity 
conservation (see I Vishniac e t al. 2003). 

Interestingly enough, while earlier on we criticized 
the concept of turbulent diffusivity as ill founded, we 
should mention that the robust reconnection of turbulent 
fields can allow parcels of magnetized fluid to move be- 
ing less constrained by the surrounding magnetic fields. 
This allows matter and heat to transfer with the rates 
determined by a simple turbulent diffusion coefhcient 
~ Vturbhurb which for typical parameters of plasma in 
galaxy cl usters exceed th e rates of thermal diffusion of 
electrons iLazarianl (200(1 ) . An important distinction re- 
mains however. Reconnection happens in thin current 
sheets, involving a negligible fraction of the magnetic 
field. Unlike a large effective resistivity, fast reconnec- 
tion cannot violate magnetic helicity conservation. 

Finally, LV99 showed that fast reconnection of stochas- 
tic magnetic field makes the models of strong MHD tur- 
bulence self-consistent. Indeed, critical balance in the 
GS95 model requires the existence of eddy-type motions 
perpendicular to the magnetic field. In the absence of 
reconnection this would result in unresolved knots that 
should drain energy from the cascade. The estimates in 
LV99 showed that the rates of reconnection predicted by 
the model are sufficient to resolve magnetic knots within 
one period. 

6.7. Final Remarks 

Another way of understanding the acceleration of energetic 
particles in the reconnection process above is to take into account 
that the length of magnetic field lines is decreasing during recon- 
nection. As a result, the physical volume of the energetic particles 
entrained on the field lines is shrinking. Thus, due to Louiville the- 
orem, their momentum should increase to preserve the constancy 
of the phase volume. 



In this paper we have considered reconnection for weak 
MHD turbulence. What would happen if the turbulence 
is strong, i.e. when the injection of turbulent energy hap- 
pens with superAlfvenic velocities? In this situation we 
expect the turbulent cascade to proceed initially in a way 
similar to ordinary hydrodynamics with strong hydrody- 
namic motions easily bending magnetic field lines. How- 
ever, at some small scale we expect the magnetic field 
to be dynamically important and resist bending. Our 
reconnection results should be applicable to turbulence 
starting at this scale. 

The successful numerical testing of the turbulent re- 
connection model presented in this paper appeals for 
more studies in this direction. It is important to stress 
that, unlike brute force simulations, we were able to test 
the LV99 model of three dimensional reconnection giving 
support to the idea that we could reliably extend our re- 
sults to much higher Lundquist numbers. The problems 
for future studies cover the importance of compressibility 
in this type of reconnection, since a substantialfraction of 
observed interstellar medium is supersonic. Another in- 
teresting question is the ability of turbulent reconnection 
to self-sustaining or to generate turbulence by itself^^. 
These are problems we hope to address in the future. 

7. SUMMARY 

In this article we put to the test the LV99 model of 
reconnection by investigating the influence of weak tur- 
bulence on the reconnection process using 3D numerical 
simulations. The turbulence is weak in the sense that 
the injection velocities are less than the Alfven speed 
and the magnetic field direction is only weakly perturbed 
by turbulence. It is strong in the sense of producing a 
strongly nonlinear cascade of energy. In our study we ex- 
perimented with different boundary conditions and with 
different ways of measuring the reconnection speed. We 
analyzed the dependence of the reconnection speed on 
the turbulence injection power, on the injection scale, as 
well as on Ohmic and anomalous resistivity. We found 
that: 

• Turbulence in 3D drastically changes the topology 
of the magnetic fleld near the interface of oppo- 
sitely directed magnetic field lines. These changes 
include the fragmentation of the current sheet, fa- 
voring multiple simultaneous reconnection events, 
as well as a substantial increase in the thickness of 
the outflow of reconnected magnetic flux and mat- 
ter. 

• The intuitive measure of reconnection defined as 
the inflow of magnetic field velocity corresponds 
well, in the condition of stationary reconnection, 
to a more rigorously defined reconnection measure 
that we introduced. The two measures show the 
same dependencies of Vrec on the power Pinj and 
injection scale Imj of the turbulence, regarding the 
different scheme used to solve MHD equations and 
different sets of initial parameters. 

^■^ While, as we argued earlier, turbulence is a natural state of 
most of the astrophysical fluids, it is important to know whether 
the reconnection is determined by the preexisting level of turbu- 
lence, or it can get self-accelerated. 
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• The reconnection rate is determined by the thick- 
ness of the outflow region. For large scale turbu- 
lence, the reconnection rate depends on the ampli- 
tude of fluctuations and injection scale as Vrec ~ 
pVj which corresponds to LV99 predictions. 

• The reconnection rate grows with the turbulence 
injection scale, which qualitatively corresponds to 
the LV99 predictions. However, in some cases the 
rate of growth is better approximated by Vrec ~ 



^inj scaling rather than Vrec ■ ~ "mj 
LV99. The difference may stem from details of the 
initial weak cascade of energy, or from limitations 
in the dynamic range available for study. 

• Reconnection in the presence of weak turbulence 
is not sensitive to Ohmic resistivity, which corre- 
sponds to the LV99 prediction that the reconnec- 
tion of weakly stochastic magnetic fields are gener- 
ically fast. 

• The introduction of anomalous resistivity does not 
change the rate of reconnection of a weakly stochas- 
tic field either, which supports the assertion in 
LV99 that in the presence of magnetic field stochas- 
ticity, which is generic for most of astrophysical en- 
vironments, plasma effects, e.g. coUisionless effects. 



lin-j predicted in 



are irrelevant in determining reconnection speeds. 

• The strength of the guide field does not change re- 
connection rates for similar rates of turbulent en- 
ergy injection. Thus fast reconnection is possible 
for generic configurations when magnetic bundles 
intersect each other at arbitrary angles. 
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